# Full text of "Number Theory Algorithms"

## See other formats

```CHAPTER 1
INTRODUCTION

This work gives a tour of the various algorithms in use for the study of the following
three important problems.

1) Given a number, determine whether it is prime or composite.

2) Given a number, factor the number completely into its prime factors.

3) Given g, t, p, find a value of x such that g x = t (mod p)

We will first consider first of these three problems.

1.1. DETERMINING WHETHER A GIVEN NUMBER IS PRIME

It is a celebrated theorem of Euler that if GCD(a, p) =1, then a <p< - p - ) =1 (mod p),
where (p is the Euler' s Phi Function, where (p(p) is the number of numbers less than p that
are relatively prime to p. When restricted to primes, (p(p) = p-1 and Euler's theorem
becomes the celebrated Fermat's Little Theorem, which states a p_I =1 (mod p), when a is
also relatively prime to p.

To get immediately an idea of whether a given number is prime or composite, we
can use the Fermat's Little Theorem. It will give the output in polynomial time, but it is
only a probabilistic algorithm. Many composite numbers will pass the test, and they are
called pseudo-primes.

In fact, to get more accurate results, more than one base can be tested, but
Carmichael Numbers pass the test for any base. An improvement to the Fermat's Little
Theorem is the Rabin Miller's Test, in which only fewer composites pass it.

1.2 PRIMALITY PROVING

In order to ensure whether the given number is prime, we have to prove it.
Unfortunately, till today there exist no general algorithm to prove easily whether the
given number is prime or not. But, for special numbers, there exist classical tests.

1 97

In 1876, Lucas used Lucas Sequences to prove that 2 -1 is prime, and this was
the largest prime number that stood for 75 years, until in 1951, Miller proved that
(2 148 +1)/17 is prime, by using Proth's Theorem. Now this test for Mersenne Numbers has
been enhanced by Lehmer, and it is called Lucas Lehmer Test. With the advent of
electronic computers over the decades, new record size primes could be discovered, and
the largest prime known to date is 2 43112609 -1, found out by GIMPS in August 2008.

Similarly, Pepin's Test is used to prove primality of Fermat Numbers. Till date,

,27

the largest Fermat Number whose cofactor has been proved to be composite is 2 +1. It
was realized by Gauss in 1796 at the age of 19 that a regular polygon of N sides can be
constructed by using straightedge and compass if and only ifN is a power of 2 times 4 or
product of one or more distinct Fermat primes. Proth's Theorem can be used for
Sierpinski Numbers and Lucas Lehmer Riesel for Riesel Numbers.

For general numbers, there are no easy tests, but there are the hard ones. The main
such tests are APRT-CLE and Elliptic Curve Primality Proving.

In August 2002, two IIT Kanpur students Kayal and Saxena, and their professor
Manindra Agrawal introduced a deterministic polynomial time algorithm for primality
proving, that does not lie on any unproved assumptions such as the Riemann Hypothesis.
[14] This algorithm is of great theoretical value, since it confirms that PRIMES is in
complexity class P, but is of questionable practical value.

Indeed, there do exist polynomial time tests for proving primality for general
numbers based upon unproved assumptions. The Riemann Hypothesis is an important
unproved conjecture in Analytical Number Theory. It states that all the non-trivial zeros
of the Zeta Function in the critical strip < Re(z) < 1 lie precisely on the line Re(z) = Vi.
If proven, many bounds on prime estimates can be improved and primality proving can
be simplified.

Assuming that the Extended Riemann Hypothesis is true, if a number N passes the
Rabin Miller's Test for all bases from 1 to 2 (log2 N) 2 , then N is prime. The Clay Math
Institute offers a \$1000000 prize for a proof of the Riemann Hypothesis.

1.3 INTEGER FACTORIZATION

Turning to the second problem, namely factoring integers into their prime factors,
Trial Division is the most obvious algorithm to factor an integer where we run through all
the prime numbers below its square root. For Mersenne Numbers 2 P -1, factors will

always be of the form 2kp+l, where p is prime, and for Fermat Numbers 2 +1, factors
will always be of the form k.2 p+2 +l. But, for general numbers, indeed, there do exist
better algorithms.

In 1970, it was barely possible to factor hard 20 digit numbers. Then, in 1974,
John Pollard introduced two new algorithms: Pollard's Rho and p-1. These algorithms are
used to find small factors, upto 15 digits in size, and are based upon the factor, not the
input.

General purpose factoring algorithms are based upon the input. In 1980, Morrison
and Brillhart introduced a sub-exponential factoring algorithm, the Continued Fraction
Method, by which it was becoming feasible to factor 50 digit numbers. [6] In 1990, Carl
Pomerance introduced the Quadratic Sieve, due to which it became feasible to factor 100
digit numbers, with the record being 1 16 digits.

In 1994, a big team comprising of over 800 people, all over the Internet, factored
the 129 digit RSA challenge number, by using the Quadratic Sieve, which had been
estimated in Martin Gardner's 1976 Scientific American Column to be safe for 40
quadrillion years. But, now Quadratic Sieve is no longer the champion. It has been
replaced by the Number Field Sieve in the spring of 1996, when it was able to factor
RSA-130 in only 15% of the time the Quadratic Sieve would have taken. [1]

With the power of computation, large numbers could be factored easily. The
record factorization by using the Special Number Field Sieve is 2 1039 -1 by Kleinjung,
et.al. in May 2007. [18] Similarly, for algorithms based upon the factor, the Elliptic
Curve Method is the algorithm in use, and it can find out factors upto 70 digits in size.
The record is a 73 digit factor from 2 -1, found out by Joppe Bos, Thorsten Kleinjung,
Arjen Lenstra and Peter Montgomery in March 2010.

In 1994, Peter Shor introduced a factoring algorithm for quantum computers,
which can factor an integer in polynomial time. In 2001, a couple of workers at IBM
factored 15 on a 4-qubit quantum computer by using this algorithm. [15]

Today, a large number of distributed computing projects are based upon the
problem of factoring integers. To name a few: Cunningham Project [31], Aliquot
Sequence [25], Home Prime Search, Fermat Numbers Factorization, ECMNET,
NFSNET, Odd Perfect Number Search, Eleven Smooth, Factoring Smarandache
Numbers, Near Repdigits, etc.

The first 11 Fermat Numbers are completely factored. In 1772, Euler showed that
641 is a factor of F5, which Fermat had previously conjectured that all such numbers are
prime. In 1880, Landry found out the factor 274177 from Fg. [11] In 1970, Morrison and
Brillhart factored F7 by using the Continued Fraction Method. [6] In 1980, Brent and
Pollard factored Fg by using a variant of Pollard's Rho. [2] Then, in 1988, Fn was the
next number to be completely factored, owing to the small size of its penultimate prime
factor. [32] In 1990, Lenstra, Lenstra, Manasse and Pollard, along with a large team
throughout the Internet, factored F9 by using the Number Field Sieve. [13] After that, F10
was the most wanted number. It was subsequently factored by Brent in 1995 by using the
Elliptic Curve Method, with a 40 digit factor. [4]

1.4 DISCRETE LOGARITHM PROBLEM

The third of the problem, called the discrete logarithm problem, comprises of
finding out a value of x, given g, t, p such that g x = t (mod p). It can be done by brute
force, or Pollard's Rho, or Pollard's Lambda or Baby Steps Giant Steps. All these have
exponential running time, and are described later. A more efficient approach to the
discrete logarithm problem is the Pohlig Hellman Algorithm, whose running time
depends upon the largest prime factor of (p(p). Index calculus method, which comprises
of combining smooth relations, has sub-exponential running time.

1.5 APPLICATIONS OF PRIME NUMBERS

Prime Numbers do have applications in cryptography, which is the art of sending
messages securely over a network. Since a message in its transit is vulnerable to attack by
intruders, a message is always encrypted and sent, and it is decrypted by the receiver. An
encryption algorithm should be secure such that it is impossible for the intruder for
decrypting the message, and read it, with his resources.

Classical algorithms for encryption and decryption are based upon the complexity
of factoring integers, with no small factors. In 1974, Diffie-Hellman key exchange was
introduced as the first public key cryptosystem. In 1978, Rivest, Shamir and Adleman
introduced a new cryptosystem and digital signature algorithm.

We choose two large prime numbers p, q. And compute their product N = pq.
(p(N) = (p-l)(q-l). Let the message be a number x < N. Each of the sender and receiver
choose their public key as an integer e, relatively prime to cp(N). They calculate their
private key by using d = e 4 (mod (p(N)). The sender encrypts the message by using
receiver's public key, with y = x e (mod N). The receiver decrypts it by using his private
key. x = y (mod N). If the intruder can factor N with his resources, then he can calculate
the private key from the public key, and he can decrypt the message. The largest RSA
key size to be factored is 768 bits, by Thorsten Kleinjung, et. al. in Bonn University,
Germany, in December 2009. This is a 232 digit number, that splits up into two 116 digit
prime factors.

Other popular crypto systems are the El Gamal cryptosystem, which is based
upon the Discrete Logarithm Problem, the problem of calculating y, given x and z such
that z = x y (mod p), and the Elliptic Curve cryptosystem. Indeed, the Elliptic Curve
cryptosystem with a smaller key size is more secure than the RSA cryptosystem of the
same key size.

Algorithms

^

J

i

1

1

Prime Testing

C

Factoring

\

Discrete Logarithms

\. J

\

j

v> J

1

1

1

Probabilistic

Deterministic

"n

r

Small Factors

>

/■

Big Factors

Brute Force

^ J

v

j

^

j

>. j

1

Fermat's Little

Special Numbers

f

General Purpose

>

(

Trial Division

\

Continued Fraction

C \

Pollard's Rho

Theorem

Method

v )

^ j

V

J

\

/

V J

v. J

Rabin Miller's Test

Lucas LehmerTest (for

APRT-CLE

>

/■

Pollard's Rho

>

f >

Pollard's Lambda

Mersenne Numbers)

v )

v J

v.

j

V

j

^ J

V )

Pepin's Test (for

f
Elliptic Curve Primality

f

Pollard's p-1

\

Number Field Sieve

Baby Steps, Giant

Fermat Numbers)

Proving

Steps

^ J

^^^^^^^^^^^B

j

\

j

^ J

V. J

Proth's Theorem (for

AKS algorithm

\

(

William's p+1

\,

Shor's algorithm

Pohlig Hellman

Sierpinski Numbers)

Algorithm

v. J

v.

j

\-

>

v J

v. J

Lucas Lehmer Riesel

Ell

ptic Curve Method

Index Calculus Method

Test (for Riesel

I ■• J

V

j

v

)

Classification of Algorithms for Integer Factorization, Primality Testing and Discrete Logarithms

CHAPTER 2
INTEGER FACTORIZATION ALGORITHMS

2.1 TRIAL DIVISION

Most obvious algorithm for integer factorization is trial division. To factor an
integer N, we test by all the prime numbers till its square root. Indeed, to speed up the
computations, we can test by 2, and then by all the subsequent odd numbers.

x = 3

f Start

)

I

Input N

Print x

Yes

Pi hit H

N = N. x

J

Print 2

I

N = N/2

No

Yes

<7 S7

Figure 2.1: A simple integer factorization algorithm

Trial Division can be used to return the small factors, but is not feasible for larger
factors. For larger factors, we need to look for better algorithms.

Time Complexity: In the worst case, we have to look for factors upto Vn, so the time
complexity is 0(vN).

2.2 POLLARD'S RHO

In 1974, John Pollard introduced a new integer factorization algorithm. Pollard's
Rho algorithm is executed on the assumption that Trial Division has eliminated factors
upto 10 5 to 10 6 . However, in the current form, one run of this algorithm returns only one
factor. More factors can be extracted by using repeated runs of the algorithm. No problem
that it is not expensive to use multiple runs for returning multiple factors.

In trial division, to factor N, we test by 2, and then examine by all the subsequent
odd numbers below Vn, in order to find out a prime factor p. That is, in the group {0, 1,
2,..., p-1}, all the factors are examined sequentially, by using a linear function. If a
quadratic function is used, say f(x) = x +a (mod p), starting with some initial value s in
the group, then we need to traverse for only lesser number of values below p, and it will
catch in a cycle sooner, and the values (mod p) will repeat.

Since we don't know what p is, we do the calculations (mod N). If two values U
and V that are the same (mod p) and different (mod N), then calculating gcd(V-U, N) will
return p.

Floyd's cycle finding method: To find a cycle, we don't need to examine all the
previous values of U and V. The tortoise hare algorithm (also known as the Floyd's cycle
finding algorithm) is an excellent algorithm. First we initialize Uo = s, Vo = s and then
compare Ui with V2, U2 with V4, U3 with W^, etc. Uk with V2k- Thus, each time we
increase the number of iterations between U and V by 1. We use the relationship Uj+i =
f(Ui) and V i+ i = f(Vi).

Since computing GCD is expensive when compared to computing successive
iterations of U and V, we set aside the task of computing GCD once every 100 iterations

or so, on the assumption that trial division has already eliminated the smaller factors.
Pollard's Rho algorithm can easily return prime factors upto 10 13 or so.

The number of iterations (mod p) needed before the algorithm gets caught in a
cycle, is called the Epact of p. To make the algorithm more effective, we need to make
the epact lesser, so we need to choose the function carefully. A linear function, ax+b
(mod p), a^O, will have the Epact to be p. Even quadratic functions such as x 2 +a (mod p),
will not be efficient, say for a=0. If the function used is x -2 (mod p), and x can be

written as y+y" 1 (mod p), then the k th iteration is given by y 2 + y" 2 (mod p).

Time Complexity: The time complexity of the algorithm depends upon the length of the
Epact (mod p), as the other processes do not consume much time. In the set of p values, if
a pseudo random function is used, how many iterations do we need on the average,
before we catch in a cycle and the values will repeat? The birthday paradox states that out
of p values, we need 1.2A/p values, on the average, such that there is 50% chance that
there are atleast two values that are same (mod p). So, it is 0(\p). For a quadratic
function, like x 2 +a (mod p), the length of the Epact is conjectured to be 0(Vp In p).The
flow chart for this algorithm is given.

Alternate cycle finding algorithms: Richard Brent proposed a new cycle finding
method, that will conserve some iterations of U and V, and is estimated to be 24% faster
than the tortoise hare algorithm. [3] In the tortoise hare algorithm, we iterate U once and
V twice, each time. In Brent's cycle finding algorithm, we save U values for iterations 1,
2, 4, 8, etc., iterations numbered 2 k and for each U value stored, we compare with V for
iterations numbered from 2 k +l to 2 k+1 -l.

Such an acceleration was used to factor Fg = 2 +1, in 1980, by Brent and Pollard. [2]

Since any factor of the Fermat Number Fg will always be of the form k.2 10 +l, a function

2 io
such as x +1 (mod p) will be very effective for this type of factors. On the execution of

this algorithm for Fg, the following factor was returned: 1238926361552897. You can use
the following mnemonic to remember this factor of Fg: "I am now entirely persuaded to
employ Rho method, a handy trick, on gigantic composite numbers. "

9

ll = s

I

( Start }

V = s

Input N

H>

Choose random "a"
between 1 and N-3

Define f(x) =
x 2 + a (mod N) j^i

Choose random "s"
between 1 and N-1

Count =

V=f(f(V))

O m = m.fV-U) mod N

m = 1

£

Old U = s H>

U = f(U)

Old V = s

Old V= V

Count = Count + 1

Print
GCD (V-U, N)

Ij

v-WflJ

/Stop "N

H>r

V=Old_V |,jj_j U=0 ld_U

~5~

U = f(U)

I

Old U = U

Pi hit
GCD (m, N)

£ No y/\s

I fGCBi

Yes

Figure 2.2: Pollard's Rho algorithm

2.3 POLLARD'S p-1 ALGORITHM

In 1974, the second factoring algorithm introduced by John Pollard is the p-1
algorithm. [29] Generally, this algorithm is performed after brute force division and
Pollard's Rho have eliminated small factors.

This algorithm gives a factor p, if p-1 is B-smooth, i.e. the largest prime factor of

p-1 lies below a bound B. Actually, the basic principle is that, by Fermat's Little

Theorem, we have that a p_1 = 1 (mod p). Further, a * = 1 (mod p), for any multiplier k.

So, if we choose a bound B (let R = LCM of all numbers below B) such that p-HR, then

a R = 1 (mod p). Computing gcd(a R -l, N) will return p.

10

We will perform the algorithm in two stages, with bounds Bl and B2. The
algorithm will return a factor p, if all the prime factors of p-1, except the largest lie below
Bl, and the largest lies below B2.
Stage 1:

Let Bl be the bound for stage 1. Let R be the LCM of all numbers below Bl.
R = qi bl .q2 b2 .q3 b3 ...qk bk , where qj are all the primes below Bl, and bi is the largest power

of q ; such that q ; b i < Bl and qi bi+1 > Bl. Let p = q 1 a i.q 2 a2 .q 3 a3 ...q k ak + 1, where q^ < Bl for
all q ; .

Now, clearly p-HR, and since a p_1 = 1 (mod p) when a, p are relatively prime, so a R = 1
(mod p). So, plGCD(a R -l, N).

The algorithm now proceeds as follows: Let R = r^^. . .rk, where rj = qi bi . Let ao = a.

r R

Calculate a; = an ' (mod N) for all i = 1 to k. We observe that ak = a (mod N), and we

evaluate GCD(ak-l, N). If all the prime factors of p-1 lie below Bl, then GCD(ak-l, N)

should return p.

Stage 2:

Suppose that all the prime factors of p-1 except the largest lie below Bl and the largest

lie below B2. So, let p = s. qi ai .q2 a2 .q3 a3 ...qk ak + 1, where p is a prime such that Bl < s <

B2. So, in this case plGCD(ak S -l, N). Let Sj, j = 1 to m, be the set of primes such that Bl <
s < B2. Put 2dj = Sj+i - Sj, and tabulate it upto the largest prime gap below s m . Since the
largest prime gap below 200000 is only 86, that below 1000000 is only 114, and that

i n

below 4.444x10 is only 652, there will not be many values to tabulate. And the big
prime gaps occur very rarely. So, there is not a need to tabulate it for all the bigger
values.

Now, that we put bi = ak Sl (mod N), and define bj + i = ak 2d ibj (mod N) for all j = 2 to m.

Now, we compute G t = GCD((b t -l)(b t+ i-l)(bt + 2-l)...(bt + c-i-l), N) for t = 1, c+1, 2c+l, ...,

[B2/c]c+l.

Since bj = a k s J (mod N), so we observe that p must divide some G t .

11

2.4 WILLIAM'S p+1 ALGORITHM

Analogous to the Pollard's p-1 algorithm, we have a William's p+1 algorithm of
integer factorization. This algorithm gives a factor p, if the prime factorization of p+1 is
(Bl, B2)-smooth, i.e. all the prime factors of p+1 except the largest lie below Bl and the
largest lie below B2. So, we need a group of order p+1. For the p-1 factoring algorithm,
the order of a (mod p) is a divisor of p-1. A group of order p+1 is based upon Lucas
Sequences. [5]
Lucas Sequences:

Let a, p be the roots of the equation x 2 -Px+Q = 0. We define U n = (<x n -|3 n )/(a-P) and V n =
(a n +p n ).

When P=l and Q=-l, then U n is the Fibonacci Sequence and V n is the Lucas Sequence.
Further, we have that U = 0, Ui = 1, and U n = PU n _i - QU n . 2 . V = 2, Vi = P, V n = PV n _i

- QV n _ 2 .

Higher U n and V n values can be calculated by using the following formula:

U 2 „ = V n U n , U 2 „-i = U n 2 - QU n .i 2 , V 2n = V n 2 - 2Q n , V 2n -i = V n V n .! - PQ nl .

Further, that U n and V n are divisibility sequences, since U n IUkn and that VJV( 2k _i) n and if

p is prime then any factor of U p will always be of the form 2kp±l and any factor of V p

will always be of the form 2kp+l, when P = 1 and Q = -1. Thus, in this case, that V p is

always congruent to 1 (mod p), and that U p will always be congruent to ±1 (mod p).

If A is given by P 2 -4Q, and p is prime, then U p -(a/ p ) = (mod p), so if (A/p) = -1, then U p+ i
= (mod p), so the order is p+1. If (A/p) = 1, then U p -i = (mod p), so the order is p-1.
Unfortunately, we don't know the value of p, and the order of group is not always p+1.
So, this algorithm returns a prime factor p, if p+1 is (Bl, B2) smooth, and (A/p) = -1, or if
p-1 is (Bl, B2) smooth and (A/p) = 1. So, it is a somewhat hybrid (p+1 and p-1)
algorithm.

The basic idea is similar to the p-1 algorithm: If (A/p) = -1, then U p+ i = (mod p), further
Uk( P +i) = (mod p) for any multiplier k > 1. So, if we take a value of R that is a multiple
ofp+l,thenplGCD(U R , N).

12

Stage 1:

Let R = r^^. . .ik as above, where rj = qi bi . So, we need to calculate ao = Ui = 1, ai = U ri ,

a2 = U ri r 2 and so on. For calculating U2n we have a formula, U2 n = V n U n = U n (PU n -

2QU n -i)- And U211-1 = U n 2 - QU n _i 2 . So, from U n we can calculate U2 n , Uztn, Ug n , etc.

But, there is a major flaw in computing Ukn from U n for the other values of k. For
calculating a 3n , we can use a 2n .a n , for other values of k: we can use modular
exponentiation for computing a kn . For example, a 59n = a 32n a 16n a 8n a 2n a n . But, for the U n
sequence, we cannot compute U3 n from U n and U2 n or other values of Ukn like this ladder.
So, we cannot compute U ri r 2 from U Fl if r2 is not a power of 2, and U ri r 2 r 3 from U ri r 2 if r3

is not a power of 2, and so on.

The only way to do this computation is to first compute the value of R at first,
and then compute Ur in a single step, by using the formula for U2n and U211-1. The
disadvantage is that R can be very large, so for larger bounds Bl, it is not practical to
store the value of R in a computer. Also, that if the bound Bl is increased to a new bound
Bl', then we have to compute the value of Ur- starting from the scratch, starting from Ui
and U2, instead of continuing it from Ur. So, we have to look for a different technique.

Let us assume that P' = (<x/p)+(p/a) and Q' = (a/p).(p/<x) = 1. Let us denote
U n (P,Q) and V n (P,Q) as the sequence generated by the roots of the equation x -Px+Q = 0.
We have U 2m (P,Q) = PQ m4 U m (P',l) (mod N). So, if p|U R (P,Q) then plU 2R (P,Q) and so
PIUr(P',1). So, we lose no generality in assuming that Q = 1. We will assume throughout
thatQ= 1.

If p is odd prime and p does not divide Q, then we have that U[p-(A/ P )] m = (mod
p), and V [p .(A/p) ]m = 2Q m[HA/p)]/2 (mod p). So, if Q = 1, then V [p . (A/p)]m = 2 (mod p).
So, if pIUr(P,1) then pl(V R (P,l) - 2). So, we can go up the ladder by using the V n
sequence.

So, we need to compute gcd(VR- 2, N), where R = r^^. . .rk. We will always use Q = 1.

13

We have V nk (P,Q) = V n (V k (P,Q),Q k ), so with Q = 1, it becomes V nk (P,l) = V n (V k (P,l),l).
So, for the stage 1 of the algorithm, we use a starting value Po such that gcd(Po -4, N) =1
We calculate P l = V ri (P ,l), P 2 = V rir2 (P ,l) = V r2 (V ri (P ,l),l) = Vr 2 (Pi,l),
P 3 = V rir2r3 (P ,l) = Vr 3 (V rir2 (P ,l),l) = V r3 (P 2 ,l) and similarly we calculate P m = V rm (P m .
i,l) (mod N) for all m = 1 to k. We observe that P k = V R (P ,1).

To calculate V r (P,l), we use the formula: V = 2, Vi = P, V 2r -i = V r V r _! - P (mod N), V 2r
= V r - 2 (mod N). We then compute gcd(P k - 2, N) to check if stage 1 returns any factor.

Stage 2:

Suppose that if (A/p) = -1 and the largest prime factor of p+1 be s. Let s be a prime such
that Bl< s < B2. Let all other factors of p+1 lie below Bl. Then, pl(V s (P k ,l) - 2) and also
plU s (P k ,l).

We tabulate U 2 d-i, U 2 d., U 2 d+i for all the distinct d, till the largest gap of successive

primes between Bl and B2, by using Uo = 0, Ui = 1, U n+ i = P k U n - U n -i.

Further, we have that U nk (P,Q) = U n (V k (P,Q),Q k ).U k (P,Q), and U m+n = U m U n+ , - QU m .

iU„, AU m+n = V m V n+ i - QV m .iV n .

So, we put that T[ Sl ] = AU s .(P k ,l) = AU s . r (Po,1)/Ur(P ,1) (mod N).

Substituting this in the above equation, we can get that

For si, we can use T[si] = P k V Sl (P k ,l) - 2V Sr i(Pk,l) (mod N) and T[s r l] = 2V Sl (P k ,l) -

PkV Sr i(P k ,l) (mod N). For other Sj, 2 < j < m, we can use T[sj+i] = T[Sj]U 2 d+i(P k -l) -

T[Sj-l]U 2dj (P k -l) (mod N), and T[s j+1 -1] = T[s j ]U 2dj (P k -l) - T[ Sj -l]U 2dj .i(Pk-l) (mod N).

We compute H t = GCD(T[st].T[st + i].T[s t+2 ]...T[s t+c .i], N) for all t = 1, c+1, 2c+l, ...,
[B2/c]c+l. Since plU s (P k ,l) if (A/p) = -1, so p must divide some H t .

14

Special Applications:

Let N = p.q where p and q are the prime factors of N. p+1 algorithm will find out
a factor p, in which p-1 algorithm has failed because it is degenerate. That is, the largest
prime factor of both p-1 and q-1 is the same. In this case, p-1 algorithm will return N,
neither p, nor q, and will fail. For example consider N =

676738418869182417666253239653933743825490524168317486197006720422834276
095213379928417101517466403166855582327009727764197997662033204334608729
42252399162065301

But, p+1 algorithm, if A is chosen such that (A/p) = 1 and (A/q) = -1, with the same
bound as that was used in the p-1 algorithm, since p-1 is smooth, and q-1 is smooth but
not q+1, the p+1 algorithm will return the factor p.

On the other hand, if we choose A such that (A/p) = -1 and (A/q) = 1, since p-1 and q-1
are smooth, and not p+1, the algorithm will return the factor q.

2.5 ELLIPTIC CURVE METHOD

Lenstra's elliptic curve factorization algorithm is an extension to Pollard's p-1
algorithm. The p-1 algorithm returns a factor p from N, if p-1 is (Bl, B2) smooth.
William's p+1 algorithm returns a factor p if either, p-1 or p+1 is smooth, depending
upon the coefficients of the recurrence equation. By the use of elliptic curves, we allow a
larger range of values for the group order, which if smooth, depending upon the curve,
will return the factor p. [22]

Elliptic Curves over Finite Fields

Elliptic Curves are defined by using the equation

y 2 = x 3 + Cx 2 + Ax + B (mod p)

This is the affine form of the elliptic curve, in which one point is assigned to the point at

infinity, O. This is the identity element of the group. The affine form can be derived from

the projective form of the elliptic curve yz = x +Cxz + Axz + Bz (mod p) in which

the point (x, y, z) of the projective form corresponds to the point (x/z, y/z) of the affine

form. The point at infinity of the projective form is given by (0, 1,0). Note that in the

15

projective form of the elliptic curve, (x, y, z) is a point on it if and only if (tx, ty, tz) is
also a point on it.

Consider the affine form of the elliptic curve. We can utilize the tremendous power of

elliptic curves only if we define some addition operation, that makes the set of points on

the elliptic curve to be an Abelian Group. We define the addition operation as follows:

Let Pi = (xi, yi) and P2 = (x 2 , y 2 ) and O be the point at infinity. Then,

-0 =

-Pi = (xi, -yi)

O + Pi = Pi

IfP 2 = -Pi, thenPi + P 2 =

If P 2 ± -Pi, then Pi + P 2 = (x 3 , y 3 ), where

X3 = m 2 - C - xi - x 2

-y 3 = m(x 3 - xi) + yi

Where the slope m is given by

m = (y 2 - yi)/(x 2 - xO if x 2 ^ xi

m = (3xi 2 + 2Cxi + A)/(2y0 if x 2 = xi

Group Order of an Elliptic Curve

We define the group order #E of the elliptic curve to be the number of points in the
elliptic curve. Since the set of points in the elliptic curve form an Abelian Group under
addition, we have [#E]P = O, where [k]P means that the point P is multiplied k times.
The factor p is returned if #E is (Bl, B2) smooth. The group order of the elliptic curve
lies between p+l-2vp and p+l+2vp.

It is the Z coordinate that determines the point at infinity, since for the point at
infinity, the Z coordinate is 0. We succeed if [#E]P = O (mod p) and not [#E]P = O (mod
N). We then calculate p = gcd(Z, N).
So, we use a form of the elliptic curve that makes use of the Z coordinate.

16

Montgomery form

The above Weierstrass form of the Elliptic Curve is not the best in terms of computations.

The Montgomery form is used instead.

gy = x + Cx + x (mod p), where A = 1 and B = 0.

For the Montgomery form, the group order is always a multiple of 4.

The Montgomery form makes use of X and Z coordinates, (X:Z) and does not matter

Suyama discovered a parameterization for which the group order is always a multiple of

12, thereby increasing the probability that it is smooth.

Take a random value of a ^ -1, 0, 1,5

Define u = o 2 - 5, v = 4a, C = [(v-u) 3 (3u+v)/4u 3 v] - 2

Then the group order will always be a multiple of 12.

We define the starting point (X:Z) to be (u :v ) mod N.

For Montgomery coordinates, we have to perform addition operation without the
use of Y coordinates. Clearly, for any given X coordinate, there will be 2 distinct values
of Y coordinate such that their sum is p. Let Pi = (Xi, Yi, Zi) and P2 = (X2, Y2, Z2). We
define X + as the X coordinate of Pi + P2, X. as the X coordinate of Pi - P2, Z + as the Z
coordinate of Pi + P2, Z_ as the Z coordinate of Pi - P2.

Then we can calculate them as follows:

(X + , Z + ) = addh(Xi, Zi, X 2 , Z 2 , X., Z.)

X + = Z.((XiX 2 - AZ1Z2) 2 - 4B(XiZ 2 + X 2 Zi + CZiZ 2 )ZiZ 2 ) and

Z + = X_(XiZ 2 - X 2 Zi) 2

Note that for Suyama's parameterization, we have A=l and B=0. So, we will have
lesser number of field multiplications for point addition, and the elliptic curve addition
operation to be homogeneous.

For doubling a point, let X+ and Z+ be the X and Z coordinates of [2]Pi. Then, we
have that (X + , Z + ) = doubleh(Xi, Zi)

17

X + = (Xi 2 - AZi 2 ) 2 - 4B(2Xi + CZOZi 3 , and
Z + = 4Zi(Xi 3 + CXi 2 Zi + AX1Z1 2 + BZi 3 ).

Elliptic Curve Multiplication

For multiplying any point P, with a scalar k, to calculate [k]P, we proceed as follows:
For example, to calculate [7]P from P, we calculate [2]P = doubleh(P),
[4]P = doubleh([2]P)

The following algorithm is used to multiply any point P (X:Z) by a scalar k.

If k = 0, then return O

Ifk= 1, then return (X:Z)

If k = 2, then return doubleh(X:Z)

If k > 2, then we define as follows:

{

(U:V) = (X:Z)

(T:W) = doubleh(X:Z)

Let k be represented as a B-bit integer in binary notation.

k B -ik B -2kB-3-..kiko

for j = B-2 down to

{

if (kj = l)then

{

(T:W) = doubleh(T:W)

}

Else if (kj = 0) then

{

18

(U:V) = doubleh(U:V)

}

}

Return (U:V)

}

Elliptic Curve Factoring Algorithm

We perform the algorithm in two stages, with bounds Bl and B2. This algorithm returns
a factor p, if all the factors of #E except the largest lie below Bl and the largest lies
below B2.
Stage 1

In stage 1, for every prime p below Bl, we take the largest power of that prime below Bl
and multiply it by the resulting point, starting with the initial point (u :v ) for the prime p
= 2. Let the Q be the resulting point after all primes below Bl have been executed. At the
end of stage 1, we take the GCD of the resulting Z coordinate with N, and check if it is
not 1 or N. If so, we return the factor p, and stop. If the GCD is 1, we proceed on to stage
2, if the GCD is N, we try out with a smaller Bl value for the same curve.
Stage 2

In stage 2, starting with the initial point Q, we multiply it by the first prime after Bl and
subsequently add it by the product of Q with the difference between successive primes
between Bl and B2. We then take the GCD of the Z coordinate of the resulting point,
with N, every 1000 iterations or so, depending upon the value of D. The following
algorithm illustrates the stage 2 of the algorithm:

Si = doubleh(Qi)
S 2 = doubleh(Si)
For d from 3 to D

{

S d = addh(S d .i, Si, S d . 2 )

p d = X(Sd) * Z(S d ) mod N

19

}
g=l

B=B1-1, if Bl is even; B should always be odd

T=[B-2D]Q

R=[B]Q

For r = B to B2; r = r+2D

{

a = X(R) * Z(R) mod N

For each prime q between r+2 and r+2D

{

5 = (q-r)/2

g = g((X(R) - X(S 5 ))(Z(R) + Z(S 8 )) - a + p d ) mod N

}

(R, T) = (addh(R, S D , T), R)

}

g = gcd(g, N)

If 1 < g < N, then return g

Else, change the a value, and try some other curve.

The following table gives the number of curves and the optimal values of Bl and
B2 needed to uncover an M-digit factor p, from N.

20

Digits

Bl

B2

Curves

20

11000

1873422

77

25

50000

12746592

206

30

250000

128992510

401

35

1000000

1045563762

948

40

3000000

5706890290

2440

45

11000000

35133391030

4590

50

43000000

240490660426

7771

55

110000000

776278396540

17899

60

260000000

3178559884516

43670

65

850000000

15892628251516

69351

Table 2.1: Optimal Bl and B2, number of curves for ECM

For example, if you want to find out factors upto 20 digits long, then you can try
77 curves with Bl = 11000 and B2 = 1873422. If it fails, the probability that there
remains a factor of 20 digits is about 1/e. If you perform 154 curves, with the same Bl
and B2 values, the probability to miss a 20 digit factor is about 1/e , and so on. So, you
can increase the Bl and B2 values to Bl = 50000 and B2 = 12746592, and then thus
perform 206 curves with these Bl and B2 values, and so on.

21

CHAPTER 3
GENERAL PURPOSE FACTORING ALGORITHMS

General Purpose Factoring Algorithms are based upon Fermat's Factorization Method.

3.1 FERMAT'S FACTORIZATION METHOD

Fermat's Factorization Method is the process of factoring a number by expressing
it as the difference of two squares. Since every composite number can be expressed as the
difference of two squares, it can subsequently be factored. If a composite number is the
product of two numbers N = pq, then it can be written as (a+b)(a-b) = a 2 - b 2 , where a =
(p+q)/2 and b = (p-q)/2.

The process is to find two numbers a and b such that a 2 = b 2 (mod N). Once, we have
such a congruence, we can calculate GCD(a+b,N) and GCD(a-b,N) to get its factors.

A Contest Problem

If someone in a contest asks to factor 8051, in this case it is useless to trial divide through
all the primes starting from 2. A quick observation will be, since 8051 = 8100 - 49 = 90
- 7 2 = (90+7)(90-7) = 83x97.

By brute force division, we would traverse forward by dividing by all the primes starting
with 2. In Fermat's Method, one would start at vN, and traverse backwards until a 2 - N is
a perfect square. Since, a general number chosen at random, has more chances to have
smaller factors, Fermat's Method is essentially worthless for a random number. It is only
generally suited for smaller numbers.

Certainly there are many enhancements to the Fermat's Method. One such enhancement
is the elimination of certain residue classes in which the factors cannot lie. For example,

22

if N = 1 (mod 4), then a cannot be even. If N = 2 (mod 3), then a must be a multiple of 3.
If a 2 - N = 17, then a = ± 1, ±2, ±4, ±8 (mod 17).

3.1.1 LEHMAN'S METHOD

Another enhancement to Fermat's Method is the search for a multiplier. We want to
search for a 2 = b 2 (mod kN), and subsequently, we can calculate gcd(a+b , N) and
gcd(a-b , N). If N = 2581, then by Fermat's Method, we have to go till a = 59, starting
with VN = a = 51, we get 59 2 - 2581 = 900 = 30 2 . But, if we use a multiplier k = 3, then
kN = 7743, and we have immediately 88 2 - l 2 = 89x87. Hence, 2581 = 89x29.

The following algorithm of Lehman formalizes the search for a multiplier.

If N has a non-trivial divisor d < N 1/3 , then return d.

For 1 < k < n 1/3 {

For 2V(kn) < a < 2^(kn) + n 1/6 /(4Vk) {

If b = v(a 2 - 4kn) is an integer, return gcd(a + b, n);

}

}

Return "N is prime".

3.2 DIXON'S FACTORIZATION METHOD

It is based upon the Fermat's Factorization Method. For larger numbers, it will be
very difficult to find out a congruence such that a 2 = b 2 (mod N). So, the idea is to find
out many squares (mod N), whose product is a square. For example, let us take N = 1649.
41 2 - 1649 = 32, 42 2 - 1649 = 115, 43 2 - 1649 = 200, none of them is a perfect square.
But their product (41.43) 2 mod 1649 = 32 x 200 = 6400 form a perfect square. So we
have 114 2 = 80 2 (mod 1649). So its factors are gcd( 114-80,1649) = 17 and
gcd(l 14+80,1649) = 97. So, we have that 1649 = 17 x 97.

The technique is to choose a collection of small primes, called the Factor Base,
and then find out values of a (mod N) whose all prime factors lie within the Factor Base
primes, i.e. they are smooth to the factor base, by using random values of a. Each such

23

congruence is called a relation, and our aim is to find out many relations whose product is
a perfect square. Since values of a 2 (mod N) are more likely to be smooth near VkN, for
integral values of k, such values are tried out.

For example, let N=l 111111 and factor base limit = 20.
We need to find values of a - N, whose all prime factors lie below 20.
For example, 1066 2 - 1111111 = 3 3 * 5 * 11* 17
1067 2 - 1111111 = 2 *3 4 * 13 2
1093 2 - 1111111 = 2 *3 3 *7* 13* 17
1121 2 -1111111 = 2*3 3 *5*7 2 * 11
1134 2 -1111111 = 5* ll 2 * 17 2
1156 2 - 1111111 = 3 2 *5 2 *7* 11 * 13

Sometimes, as the value of N grows in size, it may be very difficult to get all the
relations smooth within the factor base primes. So, we also allow partial relations, which
contains 1 or more large primes that are slightly greater than the factor base primes. Once
two relations or more have been accumulated, that contains the same large prime, we can
merge them in order to get a relation which is smooth within the factor base primes.

So, to find out a subset of these smooth relations whose product is a perfect
square, we build a matrix with the relations as rows and the factor base primes as
columns. For the factor base, only those primes p, such that the Legendre Symbol (kN/p)
= 1 are sufficient. So, different values of k generate different factor bases. The matrix
entries contain only and 1, and for every prime p, whose power is odd, we associate a
and if the power is even, we associate a 1 with it.

a

2

3

5

7

11

13

17

1066

1

1

1

1

1067

1

1093

1

1

1

1

1

1121

1

1

1

1

1

1134

1

1156

1

1

1

24

So, to find a subset whose product is a perfect square, we must find a linearly
dependent subset of rows, i.e. a set of vectors which add up to the zero vector in Z2. For
finding such a subset, algorithms such as Gaussian Elimination are available, which finds
out a linear dependent set of vectors, by bringing the matrix to an upper triangular matrix.
But, Block Lanczos algorithm [16] is faster for larger and sparse matrices. Block
Weidemann is similar to Block Lanczos algorithm, but it is suitable for parallelizing. By
elementary linear algebra, there needs to be at least M+l relations if the factor base
contains M primes, in order to guarantee that there is a linear dependency.

For example, relations {1, 2, 3, 5, 6} form a linearly dependent set of vectors. So, we
have that (1066 2 )(1067 2 )(1093 2 )(1134 2 )(1156 2 ) = 2 2 * 3 12 * 5 4 * 7 2 * ll 4 * 13 4 * 17 4 (mod
1111111). So, 62280 2 =114382 2 (mod 1111111)

Once such a congruence has been established, there is 50% chance that it will give rise to
non-trivial factorization, that is not 1 and N, or any multiple thereof.

From the above dependency, the factors are
gcd(l 14382-62280, 1 1 1 1 1 1 1) = 239 and gcd(l 14382+62280, 1111111) = 4649

3.3 CONTINUED FRACTION METHOD

In the Dixon's Factoring Method, we need to find many relations, values of x 2 -
kN, that are smooth to the factor base, so that we can find out a subset whose product is a
perfect square. Taking random relations, and trial dividing each relation by the factor
base primes is time consuming, and it is difficult to get even a single relation for larger
numbers. So, we need to find out a systematic method of finding out such relations.

Take the continued fraction expansion of VkN = qo + l/(qi + l/(q2 + l/(q3 + •••)))

Let the convergents of VkN be given by A n /B n .

The convergents can be calculated by using the following recurrence.

We have A = q , B = 1, A x = q qi+l, Bi = qi

A = qiAi.i + Ai_ 2 , Bi = qiBj-i + Bi_ 2 for all i > 2.

25

Each of the pairs of values A n and B n satisfy the congruence

A n .! 2 - kNB n .! 2 = (-l) n Q n , where q n = (VkN + P n )/Q n .

So, we have A n .i 2 = (-1)" Q n (mod N), which gives a square on the left hand side, and the

right hand side is a value, which has to be expressed as a product of primes.

Further, we always have < P n < vkN and < Q n < 2vkN, so that the Q n values are much

smaller than N, and they have high probability to be smooth.

The A n _i and Q n values can be quickly calculated by using the following formula:

Initially, we take A. 2 = 0, A_ x = 1, Q.i = kN, r_i = g, P = 0, Q = 1, g = [VkN]

And then iterate by using the following set of formula. We factor every Q n value to

determine if it is smooth. Note that we do not make use of B values.

g + P n = q n Q n + r n , where < r n < Q n

A n = q n A n .i + A n . 2 (mod N)

g + P n+ i = 2g - r n

Qn+i = Qn-i + q n (r n - r n .i)

Linear Algebra

So, we build a matrix with the relations as rows and the factor base primes as
columns. Along with the matrix, we associate a history matrix, which is a square matrix
with number of rows equal to the number of relations accumulated. Initially, the history
matrix contains 1 only in its diagonal elements. While performing Gaussian Elimination,
while adding a row to another row, the history matrix is also correspondingly adjusted.
After performing Gaussian Elimination, once a zero vector is spotted out, the
corresponding row of the history vector is analyzed to pick out the subset of relations
which when multiplied together, their product is a perfect square on the right hand side.
Square Root

Once the subset of relations is obtained, the square root procedure gives the
square root of the product of Qj's which is a perfect square. Let there be s relations Qj are

26

the various relations, this procedure calculates Q', the least positive remainder of square
root of product of Qj's (mod N).

(1) 1 = 2

(2) Q' = l

(3) R = Qi

(4) X = GCD(R,Q!)

(5) Q'=XQ'(modN)

(6) R = (R/X)(Qi/X)

(7) 1 = 1+1

(8) If I < s, then goto (4)

(9) X = VR

(10) Q' =XQ' (modN)

Now, we have that the product of various Au's = Q' (mod N), which helps to establish a

congruence of squares.

Records

2 7

In 1970, Morrison and Brillhart factored the seventh Fermat number F7 = 2 +1,
by using the continued fraction method. [6] The factor is found out to be:
59649589127497217. Use the following mnemonic to remember the factor:
First exhibitor clever John Brillhart, after numerous Saturdays, "I am allowed only
primality testing on a weekday. "

The working principle of quadratic sieve is very similar to Dixon's factorization
method, in which we find many smooth values of x - kN, which we merge together in
order to get their product to be perfect square. The working steps are the same, first we
find many smooth relations, and then solve the matrix to get the linear dependency, and
then find the square root. The last two steps are precisely the same, there is no difference

27

between them. The Dixon's Method, Continued Fraction Method, and the Quadratic
Sieve vary only in the way how the relations are accumulated.

In the Dixon's Factorization Method, we choose random values of x, and check if
they are smooth. In this way, it is very difficult to get any smooth relation as N grows
larger. So, in the Continued Fraction Method, we designed a systematic method of
accumulating relations, which come from the continued fraction expansion of vkN. For
any convergent A n /B n , we have A n 2 - kNB n 2 = (-l) n Q n , so we have that A n 2 = (-l) 2 Q n ,
which accumulates as a relation if Q n is smooth. Further, since < Q n < 2\kN, there is a
greater chance for the Q n values to be smooth.

However, in the continued fraction method, we have to trial divide every value of
Q n by the factor base primes, to check if it is smooth. This is a very time consuming
process, but we have better methods. In Quadratic Sieve, we can use a sieve to determine
those values of x, for which any factor base prime divides. The values of such x, will
follow a linear pattern, so it is very fast.

Suppose that a prime p, in the Factor Base, divides x - kN, so we will have x -
kN = (mod p), or x 2 = kN (mod p). So, the values of x for which the prime p divides x 2
- kN are given by ±VkN (mod p). This is a modular square root (mod p), and so it will
have either two solutions or no solutions (mod p), corresponding to which kN is a
quadratic residue or a quadratic non-residue (mod p). However, if kN is divisible by p,
then the two solutions x 2 = (mod p) are coincident. Note that if there are two solutions
(mod p), they are complement of each other, that is, they add up to p.

In case there are no solutions (mod p), because kN is a quadratic non-residue
(mod p), then the prime p divides no value of x - kN. So, we need not have this prime p
in the factor base. So, for choosing the factor base primes, only those primes p, for which
the Legendre Symbol (kN/p) = 1, are sufficient.

However, in case there are two solutions (mod p), they are complements of each
other, and they are given by x = ±VkN (mod p). These are the two residue classes of x for
which the prime p divides x - kN. We can sieve by every prime in the factor base, in
these two residue classes, in which two adjacent values of x, in the same residue class

28

(mod p), have a common difference of p. The advantage of sieve is that we can only test
at regular intervals for which a prime p divides x - kN, instead of testing every value of
x, if the prime p divides x 2 - kN. However, if kN = (mod p), then the two square roots
are coincident, and x = (mod p) is the only residue class, for which we have to sieve for
the values of x, for the prime p.

Within these steps, we have various enhancements. Instead of starting at x = VkN
and progressing ahead, we could have various values of x, centered at VkN. However, it
allows the values of x - kN, which are negative, but it can be managed by including -1 in
the factor base, only one more additional element in the factor base. In case the value of
x - kN is negative, then the matrix entry corresponding to -1 and this relation is 1, since
the power is assumed to be odd, and if the value of x - kN is positive, then the matrix
entry is 0, since the power is even. While merging the rows of the matrix while solving it,
to get a zero vector, the sum corresponding to -1 element should also be 0, or even, so
that the product of relations is a perfect square, which is positive.

Also that, in order for the execution of the algorithm to be fast, for every prime p

9 9

dividing x - kN, we need not store the value of x - kN, and keep on dividing it by every
prime p, as it divides x - kN. What we do is that, for every block of values of x, we use a
floating point variable for every value of x, starting from 0, and at regular intervals of p,
of the two residue classes, which we sieve over, as the prime p divides the value of x -
kN, we add the logarithm of the prime p, to the floating point variable. And at the end, if
the factor base bound is FB and the large prime bound is UB, then the smooth relations
are given by those values of x, whose floating point variables, have accumulated more
than the value of log(x 2 - kN) - log(UB).

Hensel's Lifting can be used to sieve by the prime powers. Suppose that r is a root
of the equation f(x) = (mod p m4 ), and that f(r) = (mod p m_1 ). Then we can calculate
the roots of f(x) mod p m , by using the Hensel's Lemma. Let t = [-(f(r))/p m_1 ].(f (r))" 1 (mod
p). Then, we have f(r + tp m_1 ) = (mod p m ). If x = q is a root of f(x) = (mod p m ), then x
= p m -q is the other root of this equation. Note that if f(r) = (mod p m " ) and f (r) ^ (mod
p), then f(r) = (mod p m ) has a unique solution as given above. However if f(r) = (mod

29

p m_1 ) and f (r) = (mod p), and in addition f(r) = (mod p m ), then f(r + tp m4 ) = (mod
p m ) for all integers t. But, if f(r) = (mod p m_1 ) and f (r) = (mod p), and f(r) ± (mod
p m ), then f(x) = (mod p m ) has no solution for any x = r (mod p m_I ).

Note that we use the polynomial x 2 - kN, and so we get t = [(kN - r )/p m " ].(2r)~
(mod p). But for the prime p = 2, we have to sieve separately. Since any odd square is
congruent to 1 (mod 8), so if kN = 3 (mod 4), then only 2 divides all odd values of x for
x 2 - kN. If kN = 5 (mod 8), then only 4 divides all odd values of x for x 2 - kN. If kN = 1
(mod 8), then 8 divides all odd values of x for x - kN. Sieving by higher powers of 2 can
be done by using the fact that there are 2 values of x less than 2 n_I such that their squares
are equivalent (mod 2 n ), n > 4. Their sum will be equal to 2 n4 and so we can sieve over
these values of x in the block, for 2 n , with regular intervals of 2 n " .

Suppose that kN is even, then any even square is congruent to or 4 (mod 8). So,
if kN = 2 (mod 4), then only 2 divides all even values of x for x 2 - kN. If kN = (mod 4),
then 4 divides all even values of x for x - kN. Sieving by higher powers of 2 can be done
by using the fact that the period with which we have to sieve by 2 n , n > 3 is given by
2 [(n+ ' \ where [(n+l)/2] denotes the greatest integer < (n+l)/2. We will find all the roots
of kN, that are below 2 [(n+1)/2] and then sieve over these values of x in the block for 2 n ,
with regular intervals of 2 n " . This technique can be used to sieve by the powers of 2, and
for all the powers of other primes, we can find the square root of the prime and then lift it
by using the Hensel's Lemma for the prime powers.

In this version of quadratic sieve, we use the single polynomial x - kN (mod p).
The problem with using this single polynomial is that the values of x - kN grow rapidly,
as we move away from x = "VkN. Actually, we assume that we are using the polynomial
(x+b) 2 - kN, starting with b = [VkN]+l, for positive and negative values of x. As x grows
larger, (x+b) 2 - kN ~ x 2 + 2x\kN. As x grows larger, it is less likely for the relations to
be smooth, and therefore it is more difficult to collect any relations.

30

The problem with (x+b) - kN less likely to be smooth, as x grows larger, can be
overcome by using multiple polynomials, by keeping the relations as small as possible.
Let us for simplicity purposes, keep k = 1, and we sieve polynomials of the form (ax+b)

9 9 9

- N = a x + 2axb + b - N, such that < b < a. We sieve this polynomial for x in the

9 9 9

closed interval [-M, M], Then the largest values of (ax+b) - N are given by a M - N
and the smallest values will be -N. We choose a ~ (\2N)/M to make sure that the largest
and smallest values are nearly of the same size in magnitude.

9 9 9

We choose b such that b - N is divisible by a, say b - N = ac, then (ax+b) - N =

9 9 9 9

a x + 2axb + ac, so we have that (ax+b) - N = a(ax +bx+c). Thus, every value of

9 9

(ax+b) - N will be divisible by a. Suppose that a = q , for some integer q. Then, we have

9 9 9 9

that (ax+b) = q (ax +bx+c) mod N. So, (ax+b) mod N is smooth if and only if

9 19 9

(ax +bx+c) is smooth. We have that, [(ax+b)q~ ] = (ax +bx+c) mod N. We sieve this
polynomial for all the factor base primes p, such that N is a quadratic residue (mod p),
starting with some prime q ~ V(V(2N)/M) such that N is a quadratic residue (mod q).
Then, if a prime p divides (ax+b) 2 - N, we have x = a _1 (±\N - b) mod p, and we can use
the above Hensel's Lemma to lift the roots of the prime powers. We can add the
logarithm of the prime or prime power, starting from 0, to the sieve location, if p n divides

9 9 9 9

(ax+b) - N, and since (ax+b) - N = q (ax +bx+c), so that to flag the smooth relations,
we should check for those sieve locations that have accumulated a value of atleast
log(q 2 (ax 2 +bx+c)).

The following table gives an optimal value of factor base size and value of M
needed for numbers upto 66 digits in size. [10]

Digits

Factor Base Size

M

24

100

5000

30

200

25000

36

400

25000

42

900

50000

48

1200

100000

54

2000

250000

31

60

3000

350000

66

4500

500000

Table 3.1: Optimal Factor Base Size and Block Size for MPQS

The self initiating quadratic sieve provides a fast way to change polynomials. This
uses a smaller M value than multiple polynomial quadratic sieve, and since M values are
smaller, the residues are smaller and are more likely to be smooth. This is why SIQS can
be expected to be more significantly faster than MPQS. [12]

9 9

SIQS uses the same polynomials as MPQS, namely (ax+b) - N, with b - N = ac.

9 9

So, again we have that (ax+b) - N = a(ax +bx+c). In MPQS, we took a to be a square of
a large prime. In SIQS, we take a to be the product of several factor base primes. So,

9 9

again we have that (ax+b) - N is smooth if and only if (ax +bx+c) is smooth. This
allows us to get several polynomials for each value of a, since there can be several
choices of b for every choice of a. Suppose that a = qiq2q3. . .q s , where the qi's are distinct
primes in the factor base. To generate a polynomial, we have to find out a value of b,
such that b 2 = N (mod a). There are actually 2 s distinct values of b that satisfy this
condition, but only half of the values of b will be used since (ax+b) - N gives the same
residues as (ax-b) 2 - N. So, we can use 2 s " 1 distinct polynomials for every choice of a.

We need to construct integers Bi, 1 < i < s, such that Bi 2 = N (mod qO, and that B;
= (mod qj), for 1 < j < s, i 4- j. By the Chinese Remainder Theorem, the square of any
combination of ±Bi±B2±B3...±B s will be congruent to N (mod a). Since, we want only
one value of b or -b pair, we can fix the sign of B s to be positive. Now, that the values of
B; are given by using the formula (a/qi)y aiq ., where y a , q . is the smallest positive residue of

tq^a/q;)" 1 mod qj and t q . = N (mod qi). Actually, there are two such choices for t q .. We

will choose the one that makes the value of y aiq . to be smaller.

We take the value of bi to be equal to Bi+B2+B3+...+B s (mod a). The remaining
values of b can be calculated quickly by applying a formula for Gray Code.

bi+i = b; + 2(-l) ri/2Vl B v , where 2 V divides 2i for i = 1 to 2 s4 - 1.

32

CHAPTER 4
PRIMALITY TESTS

4.1 FERMAT'S LITTLE THEOREM

The most basic method to test whether a number is prime is the Fermat's Little
Theorem. If p is prime, then a pl = 1 (mod p). We will assume throughout that a and p are
relatively prime numbers. The converse of this theorem can be used as a primality test. If
a p_1 = 1 (mod p), then p may be prime. The fact that p is prime, will not always be the
case. Though if a p_1 ^ 1 (mod p), then p is certainly composite, but if a p4 = 1 (mod p),
then p is not always prime. Composite numbers p, such that a p4 = 1 (mod p), are called
pseudo primes base a. If p satisfies a pl = 1 (mod p), then we will call p as a probable
prime (PRP) base a.

For every base a, there will be infinitely many pseudo primes, each base
containing a different set of values. For example, 341 is a PRP base 2, 91 is a PRP base 3,
217 is a PRP base 5, and 25 is a PRP base 7. So, to increase our suspicion that the
number is prime, we can test for more than one base. It would be nice if it would confirm
the number is prime.

Unfortunately, there exist numbers that are pseudo primes for every base. Such
numbers are called Carmichael Numbers. Recently it has been proved that there are
infinitely many Carmichael Numbers. Below 10 15 , there are 105212 Carmichael
Numbers, and the density increases as we go higher and higher.

But, there are Carmichael numbers with no small factors. A Carmichael number
can be recognized easily by its prime factorization. An integer N > 1 is a Carmichael
number if and only if it is composite, square free and for every prime p dividing N, we
have p-1 dividing N-l. The Carmichael Numbers below 100000 are given by 561, 1105,
1729, 2465, 2821, 6601, 8911, 10585, 15841, 29341, 41041, 46657, 52633, 62745,
63973,75361.

33

4.2 RABIN MILLER'S TEST

A better idea to check whether a number is prime, besides checking more than

one base, is to check the values of a (pl)/2 , a (p " 1)/4 , a (p " 1)/8 , etc. Note that if a (pl)/2r = 1 (mod

p), then since p is prime, a (p " 1)/2 = +1 (mod p). Otherwise, p is composite. The
advantage of this test is that fewer composite numbers pass this test.

For example, 2 170 = 1 (mod 341) and 2 85 = 32 (mod 341). So, 341 is not prime. Indeed,
we can derive its factorization from these equations. Note that 32 = 1 (mod 341), so
factors are GCD(32-1, 341) = 31 and GCD(32+1, 341) = 11. Alternately, 2 280 = 1 (mod
561) and 2 140 = 67 (mod 561). So, 67 2 = 1 (mod 561). Thus, the factors are GCD(67-1,
561) = 33 and GCD(67+1, 561) = 17.

So, we call a number p to be a strong probable prime (SPRP) base a (assuming that a, p
are relatively prime), if

(1) a*" 1 = 1 (mod p)

(2) Write p-1 = 2 s t, where t is odd. p is SPRP base a, if either a 1 = 1 (mod p) or a 2 ' = -
1 (mod p), for some i with < i < s-1 .

But, this test will not prove that the number is prime. For example, 2047 is a SPRP base
2, 121 is a SPRP base 3, 781 is a SPRP base 5, and 91 is a SPRP base 10. Every
Mersenne Number and every Fermat Number is a SPRP base 2. Although, using this test
for more than one base will give us a fairly good idea that the number is prime. But, that
is not enough. We have a few results.

If n < 1,373,653 is a both 2 and 3-SPRP, then n is prime

If n < 25,326,001 is a 2, 3 and 5-SPRP, then n is prime

If n < 118,670,087,467 is a 2, 3, 5 and 7-SPRP, then either n = 3,215,031,751 or n is prime

If n < 2,152,302,898,747 is a 2, 3, 5, 7 and 11-SPRP, then n is prime

If n < 3,474,749,660,383 is a 2, 3, 5, 7, 11 and 13-SPRP, then n is prime

If n < 341,550,071,728,321 is a 2, 3, 5, 7, 11, 13 and 17-SPRP, then n is prime

34

If n < 9,080,191 is a both 31 and 73-SPRP, then n is prime.
If n < 4,759,123,141 is a 2, 7 and 61-SPRP, then n is prime.

If n < 170,584,961 is a 350 and 395828 1543-SPRP, then n is prime.

If n < 75,792,980,677 is a 2, 379215, and 457083754-SPRP, then n is prime.

If n < 21,652,684,502,221 is a 2, 1215, 34862, and 574237825-SPRP, then n is prime.

The following number was given by Arnault in his Ph.D. thesis is a SPRP to 46 bases,

and is not prime.

803837457453639491257079614341942108138837688287558145837488917522297427

376533365218650233616396004545791504202360320876656996676098728404396540

823292873879185086916685732826776177102938969773947016708230428687109997

439976544144845341155872450633409279022275296229414984230688168540432645

75340183297861 1 1298960644845216191652872597534901

35

CHAPTER 5
PRIMALITY PROVING

The algorithms stated in the previous chapters, give an idea about whether the
number is prime, but does not prove it. In order to prove that the number is prime, we
need special algorithms, but for special numbers, polynomial time algorithms do exist.

If a p_1 = 1 (mod p), but a (p " 1)/q ^ 1 (mod p), for every prime q dividing p-1, then p is prime.
That means that the order of a (mod p) is not a proper divisor of p-1, but is equal to p-1.
Since the order of a for any number (mod p) is (p(p), this states that (p(p) cannot be less
than p-1. So, (p(p) = p-1 and that p is prime.

The problem with this to be used as a primality test is not because of finding some
value of a satisfying this, but the difficulty in obtaining the complete factorization of p-1.
As the number goes bigger, the lesser is the chance that p-1 will be smooth, and so
finding out the penultimate prime factor of p-1 will be very difficult. However, if the
factorization of p-1 is smooth, such as for Fermat Numbers in which p-1 is always a
power of 2, then this test can be used as a primality test. Note that for Fermat Numbers,
the base a can be chosen such that a is a quadratic non-residue (mod F n ), such as a = 3.
Note that only the prime 2 divides F n - 1 .

5.1 PEPIN'S TEST (FOR FERMAT NUMBERS)

F n = 2 2 " + 1 is prime if and only if 3 (F "" 1)/2 = -1 (mod F n )

The largest number tested to be composite by using this algorithm in 2000 is F24.
The above test can be extended so that we can prove the primality of a number, even only
if a partial factorization of p-1 is known.

Pocklington Theorem: If only a partial factorization of p-1 is known, such that p-1 = F.R,

where F is completely factored, and R is un-factored composite. Then if a p " = 1 (mod p)

36

and gcd(a (p ~ 1)/q - 1, p) = 1 for every prime q dividing F, then every prime factor of p is
congruent to 1 (mod F). Further, if F > Vp, then p is prime.

5.2 PROTH'S THEOREM

For numbers of the formk.2 n + 1, where factorization of p-1 is smooth, for small
values of k, the above theorem can be used. If p = k.2 n + 1 such that k < 2 n , k is odd and
a (p-i)/2 _ _ j ( moc [ p) f or some base a, then p is prime.

This theorem can be extended further such that if factorization of p-1 only upto
3 Vp is known, then it can be verified if p is prime. Let p-1 = F.R, where F is completely
factored, and R is unfactored cofactor such that Vp > F > 3 \p. Let Ci and C2 be integers
such that p = C2F 2 + CiF + 1, where < ci, C2 < F-l. Then p is prime if and only if Ci 2 -
4c2 is not a perfect square.

5.3 LUCAS LEHMER TEST (FOR MERSENNE PRIMES)

For Mersenne Numbers of the form p = 2 n - 1, one greater than it will
immediately give a smooth factorization, that is a power of 2. For Mersenne Numbers,
for which p+1 is smooth, we make use of Lucas Sequences to test if p is prime. Note that
only prime 2 divides p+1 for Mersenne Numbers.

We recall that if p is prime, then U p -(a/ P ) = (mod p), where A = P 2 - 4Q, the
recurrence equation being x 2 - Px + Q = and U n = (a n - p n )/(a - (3) and V n = (a n + p n ),
where a and P are the roots of the above equation. Similar to the Pocklington test, we can
also have a similar test. Let (A/p) = -1 and GCD(p, 2Q) = 1. If F is a divisor of p+1 and
Up+i = (mod p), and GCD(U( P+ i)/ q , p) = 1 (mod p) for every prime q dividing F, then
every prime m dividing p satisfies m = (A/m) (mod F). Further, if F > \p + 1, then p is
prime.

For Mersenne Numbers, we can use the equation x 2 - 4x + 1 = 0, and p = 2 n - 1 is

prime if and only if Y 2 n - 2 = (mod p). So, p is prime if and only if (2+V3) + (2-V3)
= (mod p). Since Vi = 4 and V2 n = V n - 2, we can use the following algorithm to test
whether p = 2 n - 1 is prime.

37

Lucas Lehmer Test: Let p = 2 n - 1, and Si = 4.
For i = 2 to n-1, define Si = S^ 2 - 2 (mod 2 n - 1).
p = 2 n — 1 is prime if and only if S n -i = 0.

5.4 LUCAS LEHMER RIESEL TEST

The Lucas Lehmer Test can be extended further to prove that numbers of the form
p = h.2 n - 1 are prime, where h < 2 n . It depends upon choosing a clever value for Si. [21]

If h = 1, then Si = 4 will hold for all odd values of n, and Si = 3 will hold for all n = 3

(mod 4)

If h = 3, then Si = 5778 can be used for all n = 0, 3 (mod 4)

If h = 6k + 1, k is an integer, then Si = (2 + Y3) h + (2 - Y3) h can be used for all the values

of n. Note that, here we use the value of P = 4 and Q = 1 and calculate the value of Vh.

Otherwise, we are left with the case where h is a multiple of 3, and we have to be more

careful in choosing the value of Si.

Try with vi values starting from 3, and let D be the square free part of vi - 4.

Let vi = a + a" 1 , and so vi satisfies the equation a 2 - via +1=0.

Now, that the values of a and r are obtained from the following equations:

a = (a + bVD) 2 /r and r = la 2 - b 2 DI, b = 1 .

Let k = V[(vi 2 - 4)/D]. Then, a = [vi + V(vi 2 - 4)]/2 = [vi + kVD]/2.

This corresponds to [(a 2 + b 2 D) + 2ab\D]/r. So, that now we can now easily solve for a

and r. These should satisfy the following equations:

(D/N) = -1, (r/N) (a 2 - b 2 D)/r = -1. If not choose a different value for vi.

Then, the starting value Si is taken to be <x h + <x" h , that is Vh with P = Vi and Q = 1,
satisfying the Lucas equation x - Px + Q = 0. The value of Vh can be calculated by using
the following recurrence equations.
V = 2, Vi = P, V n = PV„.i - QV n . 2 for all n > 2. V 2n = V n 2 - 2Q n , V 2 „-i = V n V n .i - PQ nl .

38

5.5 DISTRIBUTED COMPUTING PROJECTS

There are two distributed computing projects, one each, which make use of
Proth's Theorem and Lucas Lehmer Riesel Test for proving that numbers of the form k.2 n
+ 1 and k.2 n - 1 are prime respectively.

In 1960, Sierpinski proved that there exist infinitely many values of k for which
k.2 n + 1 is composite for all n. Such values of k are called Sierpinski Numbers. In 1962,
John Selfridge showed that k = 78557 is one such example. Seventeen or Bust is the
distributed computing project that tries to prove that 78557 is the smallest Sierpinski
Number, by finding out a value of n, for every value of k below 78557, such that k.2 n + 1
is prime. It is named so, because when the project started, there were 17 numbers below
78557, that were left to be examined. If any number exists below 78557 satisfying this
condition, then it will be the smallest Sierpinski Number.

Similarly, in 1956, Riesel proved that there are infinitely many values of k for
which k.2 n - 1 is composite for all n. These values of k are called Riesel Numbers. He
showed that k = 509203 satisfy this condition. Riesel Sieve is the distributed computing
project that tries to prove that 509203 is the smallest Riesel Number, by finding out a
value of n, for every value of k below 509203, such that k.2 n - 1 is prime. If any number
exists below 509203 satisfying this condition, then it will be the smallest Riesel Number.

39

CHAPTER 6
DISCRETE LOGARITHMS

6.1 PROBLEM STATEMENT

Given values of g, t, p, find out a value of 1 such that g 1 = t (mod p). Actually,
since g 9 ^ = 1 (mod p), so for any integral value of k, we will have g + ^ p ' = t (mod p).
So, there will be infinitely many values of 1, satisfying this condition, if there exist one,
and the problem will be to find out the lowest value of such 1. The period is the gap after
which the values oft repeat, and the period will be a divisor of (p(p). The largest value of
period will be (p(p), and the value of g, which has the period of (p(p) is called the
generator of the group (mod p), and the group is cyclic if p is prime and (p(p) = p-1. For a
generator (or primitive root) (mod p), every value oft such that 1 < t < p-1 has a solution
for 1, such that 1 < 1 < p-1.
6.1.1 Calculating the period

Calculating the period is quite obvious. First of all, you will need to factor (p(p).
Initially, assuming that the period is (p(p), then for every prime factor of (p(p), you will
need to divide that period by the prime factor. Value = period / prime factor. If you get
g Value = 1 (mod p), or g 2 * Value = g Value (mod p), if p is not prime, then the new period = old
period / prime factor. Repeat this process for all the prime factors of (p(p).

6.2 ALGORITHMS FOR SOLVING THE DISCRETE LOGARITHM PROBLEM

Assuming that on an initial attempt, any person uses brute force techniques for
small values of p, to search for a value of 1, such that g 1 = t (mod p), starting with 1=1,
we will start to look at the better algorithms to solve this problem, starting up with the
Monte Carlo Methods.

6.2.1. POLLARD'S RHO ALGORITHM

The Pollard's Rho approach can also be used to solve the Discrete Logarithm
problem. [19] Let xj = t ai g bi (mod p) and ao = 0, bo = 0, xo = 1. Here, we define a pseudo

random function such as:

40

(ai+i.bi+i) = (a ; + 1 (mod (p(p)), b ; ), if < x ; < p/3

(a i+ i, b i+ i) = (2a ; (mod (p(p)), 2b ; (mod cp(p))), if p/3 < xi < 2p/3

(ai+i, bi+0 = (a ; , b; + 1 (mod (p(p))), if 2p/3 < x ; < p

And so, we have that

xi+i = tx; (mod p), if < xj < p/3

xi+i = x; 2 (mod p), if p/3 < xj < 2p/3

xi+i = gx; (mod p), if 2p/3 < xj < p

We compare the different values of Xi to find out a match. Similar to the tortoise-
hare algorithm (Floyd's cycle finding algorithm), we compare xi with X2, X2 with X4, X3
with X6 and so on. We compare x m with X2 m each time, to look for a match. We have two
variables for x, a and b, each of which hold the iterations for m and 2m respectively.

If a match is found out, say Xj = Xk, then we have t a ig b J = t ak g bk (mod p). Since the

period is (p(p), or a divisor of (p(p), comparing the powers oft and g, and writing t = g 1 ,
we will have (aj - ak).l = bk - bj (mod (p(p)), from which we could calculate the value of 1
as to be 1 = (b k - bj)(aj - ak)" 1 (mod (p(p)).

6.2.2. POLLARD'S LAMBDA ALGORITHM

Pollard's Lambda algorithm is a distributed computing algorithm, in which
multiple computers (or processes) participate over a network. Each node starts with some
value t r for some random value of r, and uses the same pseudo random function, f: G— »S,
where S is a relatively small set of integers, whose mean value is comparable to the size
of group G. The powers g s for the values of s in set S are pre-computed. Then, the walk
starting from t r is, wo = t r , wi = wog (Wq) , W2 = Wig t(Wl) and so on. If another node starting

at r', and walking through the sequence wo', wi', W2'... has a collision with the sequence
wo, wi, W2, etc., that is Wi' = wj for some i, j, then

t-r'gfCwo'Mwj ')+••• +fl> i _ 1 ') _ j-r f(w )+f(w 1 )+...+f(Wj_ 1 )

So, since t = g , we can substitute it and calculate 1 as follows:

(r - r')l = [f(w ) + f(wi) + . . . + f( Wj .i)] - [f(w ') + f( Wl ')+... + f( Wi .i ')] (mod cp(p))

41

6.2.3 BABY STEPS, GIANT STEPS

We will look at the third algorithm for solving the discrete logarithm problem.
We have to find out a value for 1, such that g 1 = t (mod p), and whatever be the value of 1,
and if it satisfies this condition, then the lowest possible value of 1 can be given by 1 (mod
period).

Let b = Ceiling(A/p), and we will write 1 in base b = lo + lib, where < lo, li < b-
1. We have g' lb = tg"'° = th'°, where h = g" 1 (mod p). So, we will take two lists, one

containing the values of (g lb ) for < i < b-1, and the other list containing the values of
(th J ) for < j < b- 1 . We have to search for a collision, and the technique is to use some
sorting algorithm, and traverse the sorted lists to look for a collision. Sorting can be done
in worst case in 0(b In b) time by using quick sort or merge sort, and the computations of
values and collision searches can be done in 0(b) time, and so the overall time
complexity of the algorithm will be 0(b) = 0(\p). All the above three algorithms have
running time of 0(\p), and so it is exponential.

Once the collision has been detected, we have g lb = th J , so writing t = g 1 and h = g" 1 ,
we have 1 = j + ib.

6.2.4 POHLIG HELLMAN ALGORITHM

This algorithm has better running time over the above three algorithms, although
it is exponential only. The algorithm is especially fast, if the factorization of (p(p) is
smooth, and if the largest prime factor of (p(p) is, say, M, then the algorithm requires
O(wl) operations. [30] Since with the modern computers, around 10 operations are
feasible in one minute, the algorithm completes in a reasonable amount of time if the
largest prime factor of (p(p) lies below 10 1 .

Let (p(p) = qi ai .q2 a2 ...q n an , where qi, q2, ..., q n are prime numbers. We have to

calculate the value of x such that g x = t (mod p). The technique is to find out the value of
x; (mod q; a i), for i = 1,2, ..., n, and then make use of the Chinese Remainder Theorem to

find out the value of x. To find out the value of Xj (mod qj ai ), we use the following

techniques:

42

If a; = 1, and so the power of prime is 1, we calculate Xj (mod qO

We write x ; = c + Ciqi, where < Co < qi.

We get that, t^i = ( g yW% = ( g Wi%)<?(pH = ^(pHfo^v) = ( g ^Hfo ( mod p ), since

g k(p( P ) = l (mod p). We can now calculate the value of c such that t m/qi = ( s v ® , <H) c o (mod

p), in the smaller field of qi elements, by using Pollard's Rho algorithm or Baby Steps,
Giant Steps. Or we can calculate the values of (g q ») for < k < qj, if qj is smaller, and

find out the value of cq. We have that, xj = Co (mod qO

If a; = 2, the power is 2, and we need to calculate the value of Xj (mod qi )

Write X; = Co + Ciqj + C2qi 2 , where < Co, ci < qi. We can find out Co by using the above

method. We have that, t v(p)/q ? = ( g y { ^ 2 = (g c o + c i^i + W)^H 2 = (g^)V)^o.(g^) /( ii) c i

(mod p). Substituting the value of Co, we can calculate the value of ci. We observe that X;
= Co + ciq ; (mod q ; 2 ).

If a; = m > 2, then similarly calculate the value of Xj (mod q™).

Write x; = Co + ciqj + C2qi 2 + ... + c m qj m . We need to find out Co, Ci, C2, ..., c m -2 in order.

Then, we have t vW < = (g^^fo.^^p ( g ( P(P) / ii 2 ) c m-2.( g <p(p) / ii) c m-i ( mo d p) .

Substituting those already calculated values of Co, Ci, C2, ..., c m _2, we calculate the value
of Cm-i. We infer that, xj = Co + ciq; + C2qi 2 + . . . + Cm-iqi" 1 " 1 (mod q™).

Once we have the various values of Xj (mod qj ai ) for i = 1 to n, we can make use

of the Chinese Remainder Theorem to calculate the value of x. Let ri = qi ai , where we

have the various ri's for 1 < i < n, to be relatively prime. We have the different
congruences, for x = xj (mod rj). We have p = rir2r3...r n . For i = 1 to n, we calculate p; =
p/r; and v; = pf 1 (mod rj). Then, from these values, we calculate x = xipiVi + x 2 p2V 2 + ...
x n p n v n (mod p).

43

6.2.5 INDEX CALCULUS METHOD

The index calculus method for calculating discrete logarithms has sub-
exponential running time when compared to the exponential running time against the
other algorithms. The basic idea is similar to the quadratic sieve for integer factorization,
where you choose a factor base composed of small primes, and filter those values of g x
(mod p), which are smooth within the factor base primes. By elementary linear algebra, if
the factor base consists of M primes, you will need atleast M+l relations at the end, to
guarantee a linear dependency. Then, you will find a value of r, for which t r is smooth
within the factor base primes, and try to express it as a linear combination of various
values of g x . Thus, if t r = g x i +x 2 +x 3+ +x n ( mo d p ), we can calculate the value of t r = (g x ) r =
g x 1+ x 2 +x 3 +...+x N ^ mod p ^ Hence we have that, x = (xi+x 2 +x 3 +...+x N ).r" 1 (mod (p(p)).

Unlike the quadratic sieve, the difficulty with this algorithm is with the matrix. In
the quadratic sieve, we have all matrix entries to be composed of only bits containing or
1, but in the index calculus method, we have all the entries to be integers, whose size is
comparable to the value of (p(p). Clearly, we solve the matrix (mod (p(p)), or atleast (mod
Period), since the value of the period can be easily calculated earlier. So, that we have to
apply Gaussian Elimination carefully while solving the matrix and it takes extra time and
effort in solving the matrix.

44

CHAPTER 7
SUMMARY AND CONCLUSION

The algorithms Trial Division, Pollard's Rho, ECM were implemented and the time taken
for factoring numbers, with factors having upto 40 digits was studied and the results
tabulated.

Digits in Factor

Trial Division (sec)

Pollard's Rho (sec)

ECM (sec)

10

420

1

13

420000

7

15

42000000

120

1

17

3600

3

20

43200

10

25

154

30

1163

35

10296

40

72108

Table 7.1 Comparison of time taken for Trial Division, Pollard's Rho, ECM

Comparison of Trial Division, Pollard's Rho, ECM

1000000

100000

in
-a

c
o
o

CD

_cn

<D

E

10000

1000

100

10

-Trial Division
- Pollard's Rho
ECM

1 4 7 10 13 16 19 22 25 28 31 34 37 40
Digits in Factor

45

The Continued Fraction and Quadratic Sieve Methods were implemented and the time
taken for factoring numbers, with upto 55 digits was studied and the results tabulated.

Digits in Input

Continued Fraction Method (sec)

25

7

1

30

76

10

35

600

30

40

3000

95

45

21600

360

50

1800

55

10800

Table 7.2 Comparison of time taken for CFRAC with Quadratic Sieve

m
■a

c
o
o

<D

-J2-

E

10000

1000

100

10

Comparison of Continued Fraction Method with Quadratic

Sieve

Continued Fraction

1 4 7 10 13 16 19 22 25 28 31 34 37 40 43 46 49 52 55
Digits in Input

7.1 SCOPE FOR FUTURE WORK

The Quadratic Sieve is not the best known factoring algorithm for general purpose

numbers with two large prime factors, but the Number Field Sieve has better asymptotic
running time. In fact, the cross over point between the Quadratic Sieve and Number Field

46

Sieve is about 97 digits, above which the Number Field Sieve is faster than the Quadratic
Sieve. To be precise, Continued Fraction Method has running time of e (2 og n og og , the
Quadratic Sieve has running time of e ( og n og og n) and the Number Field Sieve has

running time of e c(log n) (log log n) . The Number Field Sieve uses a higher degree
polynomial, and it comes in two variations, the Special Number Field Sieve and General
Number Field Sieve. The Special Number Field Sieve can only be used for numbers of
special algebraic form, such as Cunningham Numbers [31], which are of the form b n + 1,
and in general, it should be of the form qo + qix + q2X 2 + . . . q n x n , for small values of qi, if
a degree n polynomial is used. For SNFS, degree 4 polynomial is used for 100 digits,
degree 5 for 150 digits and degree 6 is optimal for 200 digits. For example, in 1990, the

S 1 9

ninth Fermat number Fg = 2 +1 was completely factored by using the Special Number
Field Sieve. [13] If a general number can be expressed in the above form, with larger
values of q;, then it becomes the General Number Field Sieve, differing only in the
polynomial selection method, with the larger values of coefficients qj, and this algorithm
can be used for any general number. [17] In fact SNFS is about 1.43 times faster than
GNFS for a number of the same size.

Further, Shor's algorithm is a polynomial time factoring algorithm on a quantum
computer. It is based upon finding out the period p for N, by using quantum
computations, such that a p = 1 (mod N), where a is chosen such that it is relatively prime
to N. Then, we can calculate the value of (p(N), which will be a multiple of the period,
and after that we can easily extract the factors of N. [15]

For discrete logarithms, Index Calculus Method has sub-exponential running
time. For primality proving, we can work on APRT-CLE [20] and Elliptic Curve
Primality Proving. Elliptic Curve Primality Proving is based upon the fact that if s divides
k, and that [k]P = O (mod p), but [k/q]P 4- O (mod p) for every prime q dividing s, then
every prime m dividing p satisfies, [#E]F m = (mod s). Further, if s > Vp, then p is prime.
Actually, the techniques to find out the group order of the elliptic curve digs deep into
polynomial factorization, so that we will work upon it later on.

47

In August 2002, two IIT Kanpur students: Kayal and Saxena, and their professor:
Agrawal, introduced a new deterministic polynomial time algorithm for proving
primality, that does not lie on any unproved assumptions, such as the Riemann
Hypothesis. [14] It states that p is prime if and only if (x-a) p = x p - a (mod p). Here x and
a are not simple constants, but it should be treated as a polynomial variable, which is a bit
confusing, so we will analyze about it and work upon it later on.

In addition, for the already implemented algorithms, for the Linear Algebra of
Quadratic Sieve and other allied algorithms, there are better algorithms than Gaussian
Elimination for solving sparse matrices over GF(2). One of them is Block Lanczos, [16]
which can be implemented on a single computer, and Block Weidemann, which is
efficient as a distributed computing algorithm.

For algorithms such as p-1, p+1 and ECM, we have to multiply the initial point
with every prime less than Bl, and for every prime, we should use the largest power of
prime less than Bl. There are very efficient algorithms for Big Integer multiplication,
such as Fast Fourier Transform algorithms. For Big Integer squaring, such as in Lucas
Lehmer Test and Pepin's test, we can use Discrete Weighted Transform algorithms. For
these two algorithms, it is very convenient to keep 2-1 and 2 + 1, as a string of bits,
in base 2, consisting only of zeros and ones, and it is very comfortable to square or

2 n

multiply (mod 2 - 1) or (mod 2 + 1). Keeping these numbers in base 10, and doing
such calculations, in decimal would take a lot of time.

48

REFERENCES

Papers:

1. "A Tale of Two Sieves" by Carl Pomerance

2. "Factorization of the Eighth Fermat Number" by Richard Brent and John Pollard

3. "An Improved Monte Carlo Factorization Algorithm" by Richard Brent

4. "Factorization of the Tenth Fermat Number" by Richard Brent

5. "A p+1 method of Factoring" by Hugh Williams

6. "A method of Factoring and the Factorization of F7" by Michael Morrison and John Brillhart

7. "An FFT extension to the p-1 Factoring algorithm" by Peter Montgomery and Robert
Silverman

8. "Speeding the Pollard and Elliptic Curve Methods of Factorization" by Peter Montgomery

9. "A practical analysis of the Elliptic Curve Factoring Algorithm" by Robert Silverman and
Samuel Wag staff

10. "The Multiple Polynomial Quadratic Sieve" by Robert Silverman

11. "How was F6 factored?" by Hugh Williams

12. "Factoring Integers with the Self Initializing Quadratic Sieve" by Scott Contini

13. "The Factorization of the Ninth Fermat Number" by Arjen Lenstra, Hendrik Lenstra, Mark
Manas se and John Pollard

14. "Primes is in P" by Manindra Agrawal, Neeraj Kayal and Nitin Saxena

15. "Quantum Computing and Shor's algorithm" by Matthew Hayward

16. "A Block Lanczos Algorithm for finding dependencies over GF(2)" by Peter Montgomery

17. "On Polynomial Selection for the General Number Field Sieve" by Thorsten Kleinjung

18. "A kilobit Special Number Field Sieve Factorization" by Arjen Lenstra, Kazumaro Aoki,
Jens Franke, Thorsten Kleinjung and Dag Arne Osvik

19. "Monte Carlo Methods for Index Computation (mod p)" by John Pollard

20. "On distinguishing prime numbers from composite numbers" by Leonard Adleman, Carl
Pomerance and Robert Rumely

21. "Lucasian Criteria for the Primality of N = h-2 n - 1" by Hans Riesel

22. "Factoring Integers with Elliptic Curves" by Hendrik Lenstra

49

23. "The Number Field Sieve" by Steven Byrnes

24. "Modular Multiplication without Trial Division" by Peter Montgomery

25. "Aliquot Sequence 3630 ends after reaching 100 digits" by Manuel Benito, Wolfgang
Creyaufmiiller, Juan Varona and Paul Zimmermann

26. "A Beginner's Guide to the General Number Field Sieve" by Michael Case

27. "The Multiple Polynomial Quadratic Sieve" by Brandt Kurowski

28. "Factoring Large Numbers with a Quadratic Sieve" by Joseph Gerver

29. "Theorems on Factorization and Primality Testing" by John Pollard

30. "An Improved Algorithm for Computing Logarithms over GF(p) and its Cryptographic
Significance" by S.Pohlig and M.Hellman

31. "Factorizations of b n ± 1, b = 2, 3, 5, 6, 7, 10, 11, 12 up to high powers" by John Brillhart,
Dick Lehmer, John Selfridge, Bryant Tuckerman and Samuel Wagstaff, Jr.

32. "Factorization of the Eleventh Fermat Number" by Richard Brent

Books:

1. "The Art of Computer Programming" by Donald Knuth, Volume 2: Semi-Numerical
Algorithms

2. "Prime Numbers - A Computational Perspective" by Richard Crandall and Carl Pomerance

3. "The Development of The Number Field Sieve" by Arjen Lenstra, Hendrik Lenstra, Mark
Manas se and John Pollard

4. "Prime Numbers and Computer Methods for Factorization" by Hans Riesel

5. "The Arithmetic of Elliptic Curves" by Joseph Silverman and John Tate

50

APPENDIX A
NUMBER THEORY CONCEPTS

A.l GCD AND INVERSE

If x, y are two numbers, then gcd(x, y) represents the greatest common divisor of x and y.
The gcd of two numbers can be calculated by the following code. Let x, y be two numbers
such that x > y
while (y > 0)

{

(x, y) = (y, x mod y)

}
return x

Bezout's Theorem

If x, y are two numbers such that GCD (x, y) = g, then there exist integers a, b such that
ax + by = g. The constants a, b can be calculated by the following algorithm.

(a, b, g, u, v, w) = (1, 0, x, 0, 1, y)
while (w > 0)

{

q = Integer Part of (g / w)

(a, b, g, u, v, w) = (u, v, w, a - qu, b - qv, g - qw)

}

return (a, b, g)

So, if you want to calculate inverse of x (mod y), that is, the number x" 1 such that x x x" 1
= 1 (mod y), then, it is a. And the inverse of y (mod x) is b.
The inverse can be calculated more easily by using the following formula.
x" 1 (mod y) = x 9 ^ " * (mod y), where (p is the Euler's Phi Function.

51

A.2 CHINESE REMAINDER THEOREM

Consider the following set of congruences,
x = ni (mod pi)
x = n2 (mod P2)
x = n 3 (mod p 3 )

x = n r (mod p r ), where pi, P2, P3, ... p r are pairwise relatively prime, i.e. gcd(pj, pj) = 1 for

all 1 < i < j < r and LCM (p x , p 2 , p 3 , ... p r ) = piP2P3-..p r

The Chinese Remainder Theorem claims that x has a unique solution (mod pip2P3-..p r )

under these set of conditions such that < x < pip2P3-..p r

For 1 < i < r, if mj = pip2p3...p r / Pi and v; = mi" 1 (mod pO, then x is given by

x = nimiVi + n2m2V2 + n3ni3V3 + ... + n r m r v r (mod pip2P3...p r )

A.3 MODULAR SQUARE ROOT ALGORITHM

Here computing the square root of a number modulo a prime number is discussed.

Let a be the given number and p be the given prime number which forms the modulo for

the calculations. If there exists a number x such that x = a (mod p), then x is called a

square root of a (mod p)

If there exists a square root for a (mod p), then a is called quadratic residue (mod p), then

there will be two square roots for a (mod p), such that their sum will be p.

If there exist no square root for a (mod p), then a is called quadratic non-residue (mod p)

In short, half of a less than p, have two square roots, and half of a less than p, have no

square roots.

The Legendre Symbol (a/p) is defined by the following statements

1 if a is quadratic residue (mod p)

-1 if a is quadratic non-residue (mod p)

if a and p are not relatively prime.

52

For any prime number p and a base number a, the Legendre Symbol can be calculated by

using the following formula

(a/p) = a (l>1)/2 (mod p)

Note that, if p is prime number and a is less than p, by the Fermat's Little Theorem, a p_1 =

1 (mod p). Since p is prime, a (p " 1)/2 should be either 1 (mod p) or -1 (mod p).

If m is a composite number (or a prime number with only one prime factor), and has

prime factorization pi b x p 2 c x p 3 d x ... x p n z then the Jacobi Symbol (a/m) is defined by

(a/ Pl ) b x (a/p 2 ) c x (a/p 3 ) d x ... x (a/p n ) z .

So, given a base number a and a prime number p, whether a has a square root
(mod p), or not can be calculated by the Legendre Symbol (a/p) = a (pl)/2 (mod p).

If (a/p) = 1, then a has two square roots (mod p), both of which add up to p.

If (a/p) = -1, then a has no square roots (mod p)

So, if (a/p) = 1, there exist two square roots for a (mod p)

We will see and then discuss here how to compute those square roots x, given a and p.

However note that this algorithm works only when p is prime.

If p = 3 (mod 4), then it is the best case. Square root is given by the following formula

x = a (p+1)/4 (modp)

If p = 5 (mod 8), Square root is calculated with the following formula x = a (p+3)/8 (mod p)

9 9

Now, with this, either x = a (mod p) or x = (p-a) (mod p)

2

If x = (p-a) (mod p), then
Assign x = 2 (p " 1)/4 x (mod p)

If p = 1 (mod 8), then it is the worst case.

Square root can be calculated by using the following algorithm

Select a random d between and p-1 such that (d/p) = -1, i.e. d is a quadratic non-residue

(mod p). Write p-1 = 2 s t where t is odd, i.e. t is the highest odd divisor of p-1

Let A = a' mod p

Let D = d' mod p

53

Assign m =

For (i from to s-1)

{

l{((AD m f 1 ' 1 = -l(modp))

{

Assign m = m + 2 1

}
}

x = a (t + l)/2 D m/2 (modp)

return x

A.4 SUM OF FACTORS AND EULER'S PHI FUNCTION

The sum of all the positive factors of a given number (including itself) can be
found out by the given formula. If the given number m can be written as product of prime
factors as m = pi a x p 2 b x ... x p n z , then the sum of factors of m, is given by

(1 + Pl + pi 2 + ... + pi a ) (1 + p 2 + P2 2 + ••• + P2 b ) ••• (1 + Pn + Pn 2 + ••• + Pn)

= ( Pl a+1 -l)/( Pl -l) x (p 2 b+1 -l)/(p 2 -l) x ... x ( Pn z+1 -l)/(p n -l)

= (( P i a+1 -1) X (p 2 b+1 -l) X ... X (p n Z+1 -l)) / ((pi-1) X (P2-D X ... X ( Pn -1))

The Euler's Phi Function is defined as the number of numbers, below the given number,
and relatively prime to the given number. The property of Phi Function is that if a, m are
relatively prime, then a 9 ' } =1 (mod m).And if m is prime, then (p(m) = m-1, So, if p is
prime, a pl = 1 (mod p)

The Phi Function can be calculated by using the following formula.

If m is the given number m = pi a x p 2 b x ... x p n z , then (p(m) is given by

mx(l - (1/pO) x (1 - (l/p 2 ))x ... x(l - (l/p n )) = (m x ( Pl -1) x (p 2 -l) x ... x (p n -l)) / ( Pl x p 2

X ... Xp n )

54

The Mobius function of a given number is defined as 0, if it has at least one repeated
prime factor, 1, if all prime factors are distinct and it has even number of prime factors
and - 1 , if all prime factors are distinct and it has odd number of prime factors

A.5 PROPERTIES OF LEGENDRE SYMBOL AND JACOBI SYMBOL

l.(a/p) = a (p " 1)/2 (modp)

2. (ab/m) = (a/m).(b/m)

3. (a/mn) = (a/m).(a/n)

4. (-l/m) = (-l) (m - 1)/2

5. (2/m) = (-l)( m2 - 1)/8 , that is (2/m) = 1 if m = 1 or 7 (mod 8), and (2/m) = -1 if m = 3 or 5
(mod 8)

Finally, the law of quadratic reciprocity gives the value of (m/n) if (n/m) is known, and if
m and n are relatively prime. If m and n are not relatively prime, then both their value is
defined to be 0.

(m/n).(n/m) = (.i^Xn-iy^ mat is if atleast one of m, n = 1 (mod 4), then (m/n) = (n/m).
If both of m, n = 3 (mod 4), then (m/n) = -(n/m).

A.6 THE ALIQUOT SEQUENCE

Mersenne Primes are closely associated with Perfect Numbers. Every new discovery of a
Mersenne Prime generates a new Perfect Number. A perfect number is a number whose
sum of all its factors except itself is equal to the same number. The first few perfect
numbers are 6, 28, 496, 8128, 33550336, 8589869056, 137438691328,
2305843008139952128. It is not difficult to prove that if 2 P - 1 is prime, then (2 P " 1 )(2 P -1) is
a perfect number and that every even perfect number has this form. Every even perfect
numbers other than 6 will always be congruent to 1 (mod 9), and will end in 6 or 28.

It is not known whether any odd perfect numbers exist, but it has been proved that there
aren't any less than 10 300 . Carl Pomerance has given a heuristic estimate that odd perfect
numbers are extremely unlikely. If any odd perfect number exists, then it must be a
perfect square times an odd power of a single prime, it must be divisible by atleast 9

55

distinct primes and has atleast 75 prime factors (not necessarily distinct), and then it has a
prime factor greater than 10 .

Let us define the function s(n) to be the sum of all factors of n except itself, and a(n) to
be the sum of all factors of n including itself. To be more clear, s(n) = o(n) - n. For a
perfect number, we have that s(n) = n, or a(n) = 2n. This function s(n) is called the
Restricted Divisor Function. If s(n) > n, then n is called as Abundant Number, and if s(n)
< n, then n is known as Deficient Number. For example, 945 is the least Odd Abundant
Number, and that any multiple of 6 will always be abundant.

There exist plenty of pair of numbers (x, y) such that s(x) = y and s(y) = x. Such pairs of
numbers are called Amicable Numbers.

The amicable pairs below 20000 are given by: (220, 284), (1184, 1210), (2620, 2924),
(5020, 5564), (6232, 6368), (10744, 10856), (12285, 14595), (17296, 18416).

Let us iterate the Restricted Divisor Function. Let s'(x) be the value obtained after i
iterations of the Restricted Divisor Function, starting with x. By this process of iterating
the Restricted Divisor Function, we obtain a sequence. That sequence is called the
Aliquot Sequence. Suppose that s k (x) = x, i.e. after k iterations of the Restricted Divisor
Function, starting with x, we again obtain x, then that forms an Aliquot cycle of order k.
In order words, x, s(x), s 2 (x), s 3 (x), ..., s k4 (x) are called sociable numbers of order k,
where k> 3.

For example, (12496, 14288, 15472, 14536, 14264) forms an Aliquot cycle of order 5.
The largest known Aliquot cycle is of order 28:

(14316, 19116, 31704, 47616, 83328, 177792, 295488, 629072, 589786, 294896,
358336, 418904, 366556, 274924, 275444, 243760, 376736, 381028, 285778, 152990,
122410, 97946, 48976, 45946, 22976, 22744, 19916, 17716)

56

There exist plenty of Aliquot cycles of order 4, most of them can be generated by using
Thabit rules. The first ten of them being:

(1264460, 1547860, 1727636, 1305184), (2115324, 3317740, 3649556, 2797612),
(2784580, 3265940, 3707572, 3370604), (4938136, 5753864, 5504056, 5423384),
(7169104, 7538660, 8292568, 7520432), (18048976, 20100368, 18914992, 19252208),
(18656380, 20522060, 28630036, 24289964), (28158165, 29902635, 30853845, 29971755),
(46722700, 56833172, 53718220, 59090084), (81128632, 91314968, 96389032, 91401368)

We also know five Aliquot cycles of order 6:

(21548919483, 23625285957, 24825443643, 26762383557, 25958284443,

23816997477)

(90632826380, 101889891700, 127527369100, 159713440756, 129092518924,

106246338676)

(1771417411016, 1851936384424, 2118923133656, 2426887897384, 2200652585816,

2024477041144)

(3524434872392, 4483305479608, 4017343956392, 4574630214808, 4018261509992,

3890837171608)

(4773123705616, 5826394399664, 5574013457296, 5454772780208, 5363145542992,

5091331952624)

Two Aliquot cycles of order 8:

(1095447416, 1259477224, 1156962296, 1330251784, 1221976136, 1127671864,

1245926216, 1213138984)

(1276254780, 2299401444, 3071310364, 2303482780, 2629903076, 2209210588,

2223459332, 1697298124)

and one Aliquot cycle of order 9:

(805984760, 1268997640, 1803863720, 2308845400, 3059220620, 3367978564,

2525983930, 2301481286, 1611969514)

That is about Aliquot cycles. Coming back to Aliquot Sequences, let us take some
examples starting with arbitrary values.

57

12, 16, 15, 9, 4, 3

30, 42, 54, 66, 78, 90, 144, 259, 45, 33, 15, 9, 4, 3

60, 108, 172, 136, 134, 70, 74, 40, 50, 43

96, 156, 236, 184, 176, 196, 203, 37

120, 240, 504, 1056, 1968, 3240, 7650, 14112, 32571, 27333, 12161

48, 76, 64, 63, 41

24, 36, 55, 17

18,21,11

Let us now consider with a longer sequence

180, 366, 378, 582, 594, 846, 1026, 1374, 1386, 2358, 2790, 4698, 6192, 11540, 12736,

12664, 11096, 11104, 10820, 11944, 10466, 5236, 6860, 9940, 14252, 14308, 15218,

10894, 6746, 3376, 3196, 2852, 2524, 1900, 2440, 3140, 3496, 3704, 3256, 3584, 4600,

6560, 9316, 8072, 7078, 3542, 3370, 2714, 1606, 1058, 601

Another example

11025, 11946, 14262, 14274, 19578, 22758, 22770, 44622, 56154, 75174, 101082,

113190, 232410, 338982, 450354, 470094, 490674, 509838, 680562, 844764, 1314372,

1952108, 1496764, 1132100, 1324774, 843074, 428734, 228194, 119134, 59570, 71758,

35882, 31510, 28106, 20278, 10142, 6490, 6470, 5194, 4040, 5140, 5696, 5734, 3194,

1600, 2337, 1023, 513, 287, 49, 8, 7

All these sequences terminate in a prime. Some sequences can terminate in a perfect

number, or an amicable pair, or Aliquot cycle of higher order as well. For example:

2195,445,95,25,6,6...

650,652,496,496...

2856, 5784, 8736, 19488, 40992, 84000, 230496, 475356, 792484, 1013852, 1013908,

1058092, 1264340, 2049964, 2123576, 2778664, 3492536, 3077104, 2884816, 3391568,

3775384, 3303476, 3003244, 2288756, 1773484, 1558964, 1440718, 720362, 360184,

376736, 381028, 285778, 152990, 122410, 97946, 48976, 45946, 22976, 22744, 19916,

17716, 14316, 19116, 31704, 47616, 83328, 177792, 295488, 629072, 589786, 294896,

358336, 418904, 366556, 274924, 275444, 243760, 376736, ...

Now try that out with some longer sequences as well:

58

s 176 (138) = 59, s 299 (702) = 191, s 193 (720) = 277, s 183 (936) = 43, s 800 (4488) = 43,
s 671 (1848) = 43, s 177 (976950) = 6, s"(980460) = 2620, s 228 (17490) = 1264460, s 393 (2880)
= 277, s 746 (840) = 601, s 1074 (1248) = 37, s 166 (858) = 59

In 1888, Eugene Catalan conjectured that all Aliquot Sequences terminate in a prime, or
perfect number or amicable pair or Aliquot cycle of higher order. This is one of the
greatest problems of computational number theory. It is unlikely that this conjecture is
true, for there are some sequences that grow indefinitely without bound, without even
catching up in a cycle. These sequences are called open-end sequences. According to the
current factoring resources, these sequences have been extended atleast upto 120 digits in
size.

Today, there are twelve open-end sequences below 1000:

276, 306, 396, 552, 564, 660, 696, 780, 828, 888, 966, 996

But, notice that (276, 306, 396, 696) all belong to the same family. And so does (660,

828, 996), (564, 780) and (552, 888)

Considering the smallest member of distinct sequences, we have 5 open-end sequences

below 1000, namely 276, 552, 564, 660, 966. These are called as "Lehmer Five". It was

originally known as "Lehmer Six" along with 840, which had been later on terminated by

several groups separately.

The open-end sequences between 1000 and 2000 were termed as "Godwin Fourteen"

including 1248 and 1848, but today there are only twelve of them: 1074, 1134, 1464,

1476, 1488, 1512, 1560, 1578, 1632, 1734, 1920, 1992.

Today, there are exactly 81 open-end sequences below 10000.

276, 552, 564, 660, 966, 1074, 1134, 1464, 1476, 1488, 1512, 1560, 1578, 1632, 1734,
1920, 1992, 2232, 2340, 2360, 2484, 2514, 2664, 2712, 2982, 3270, 3366, 3408, 3432,
3564, 3678, 3774, 3876, 3906, 4116, 4224, 4290, 4350, 4380, 4788, 4800, 4842, 5148,
5208, 5250, 5352, 5400, 5448, 5736, 5748, 5778, 6396, 6552, 6680, 6822, 6832, 6984,
7044, 7392, 7560, 7890, 7920, 8040, 8154, 8184, 8288, 8352, 8760, 8844, 8904, 9120,
9282, 9336, 9378, 9436, 9462, 9480, 9588, 9684, 9708, 9852
Any progress in computations can certainly be able to reduce these numbers down.

59

A.7 WIEFERICH PRIMES

Any prime number p satisfies with a p_1 = 1 (mod p) for any base 'a' such that GCD(a, p)
= 1. Thus, 2 pl = 1 (mod p) for all primes p > 3. A Wieferich prime is a prime p which
satisfies with 2 P_1 = 1 (mod p 2 ). If square of some prime q divides a Mersenne number
with a prime exponent p, i.e. q 2 I 2 P -1, where p and q are primes, then q must be a
Wieferich prime. The Wieferich prime concept was first introduced by Wieferich in
1909, that if the first case of Fermat's Last Theorem is false for prime p, then p must be a
Wieferich prime. Only two Wieferich primes are known till today, 1093 and 3511, found
out by Meissner in 1913 and Beeger in 1922 respectively. This search has been carried

To answer whether there exist any prime q such that q 2 I 2 P -1, where p is also prime. 1093
does not satisfy the necessary condition, i.e. 2 (p " 1)/2 = 1 (mod p 2 ), that is also needed for
the above statement to be true. 3511 does not divide any Mersenne number with prime
exponent at all, so no such prime q is being known till today.

Similar to Wieferich primes, that satisfy with 2 P ~ 1 = 1 (mod p 2 ), that are based upon the
Fermat's Little Theorem a p_1 = 1 (mod p). There is also some research carried upon
Wilson primes, that satisfy with (p-1)! + 1=0 (mod p 2 ), that is based upon the Wilson's
Theorem (p-1)! + 1=0 (mod p). And then Wall-Sun-Sun primes, that satisfy with F p .( P /5)
= (mod p 2 ), that is based upon the Fibonacci number primality test, that is F p .( P /5) =
(mod p). Here Fk denotes the k th Fibonacci number, and that (p/5) is the Legendre
Symbol. The only known Wilson primes are 5, 13 and 563 (below 500 million), and then
no Wall-Sun-Sun prime is being known till today (below 200 trillion).

A.8 FACTORS OF SOME SPECIAL NUMBERS

Some people believe that if q = 2 P -1 is prime, then 2 q -l is also prime. However, that is
not true. But if q = 2 P -1 is composite (with some prime factor F), then 2 q -l will always be
composite (with 2 F -1 as a factor). For the first four Mersenne primes 2, 3, 5 and 7,

60

however 2-1,2-1,2 -1 and 2 -1 are all prime. But for the next four Mersenne primes

13, 17, 19 and 31, we know that 2 8191 -1, 2 131071 -1, 2 524287 -l, 2 2147483647 -1 all are composite.

2 8191 -1 has factors 338193759479, 210206826754181103207028761697008013415622289

2 131071 -1 has factors 231733529, 64296354767

2 524287 -l has factors 62914441, 5746991873407, 2106734551102073202633922471,

824271579602877114508714150039

2 2147483647 -1 has factors 295257526626031, 87054709261955177,

242557615644693265201, 178021379228511215367151

Tony Forbes is leading a project to search for a factor of the next term 2 2305843009213693951 -1,

none being found out so far.

Corning to the case of Fermat primes, 2-1,2-1,2 -1 all are primes. But the next term
2 257 -l is not prime. 2 257 -l = 535006138814359 x 1155695385246619182673033 x
3745505985018 1093658 17766300963 13 18 1393
2 65537 -l has factors 513668017883326358119, 883852340565787164089923172087

The first few Mersenne numbers, that are known to be composite, but with no known
prime factors at all, till today are 2 1061 -1, 2 1277 -1, 2 1619 -1, 2 1753 -1.

A.9 SUM OF SQUARES

Sum of two squares:

Every prime number of the form 1 (mod 4) can be uniquely represented as sum of two
squares. For example, 5 = 2 2 + l 2 , 13 = 3 2 + 2 2 , 17 = 4 2 + l 2 , 29 = 5 2 + 2 2 , 37 = 6 2 + l 2 , 41
= 5 2 + 4 2 , 53 = 7 2 + 2 2 , 61 = 6 2 + 5 2 , 73 = 8 2 + 3 2 , 89 = 8 2 + 5 2 , 97 = 9 2 + 4 2 .

9 9

2 has also unique representation: 1 + 1 . Any number of form 3 (mod 4) cannot be
represented as sum of two squares in any way. This algorithm can be used to return the

9 9

values (a, b) for a +b representation for any arbitrary prime number p of form 1 (mod 4).

Initialize a = p, choose any quadratic non-residue (mod p) and then assign its value to
variable x. Initialize b = x (p " 1)/4 (mod p).

61

While (b 2 > p)

{

(a, b) = (b, a mod b)

}

Return (b, a mod b)

If any integer N has a prime factor of form 3 (mod 4) to an odd power, then N cannot be
represented as sum of two squares. If N has all prime factors of form 3 (mod 4) to even
powers (>2), then N has no sum of two squares representation in lowest terms.
That is, if N = y 2k z, where y 2k is product of all prime factors of N of form 3 (mod 4), and
that z is product of all prime factors of N of form 1 (mod 4), then for every sum of two

9 9

squares representation of z = q + r , the only possible representations of N as sum of two
squares are N = (y k q) 2 + (y k r) 2 -

If N is the product of two distinct primes of form 1 (mod 4), then N will have two

different representations as sum of two squares. This is demonstrated by using the

following equation:

(a 2 + b 2 ) (x 2 + y 2 ) = lax + byl 2 + lay - bxl 2 = lax - byl 2 + lay + bxl 2

This can be similarly extended when N is product of multiple primes of form 1 (mod 4)

not necessarily distinct. If N is the product of k distinct primes of form 1 (mod 4), then N

can be represented as sum of two squares in 2 k_1 ways. If the prime factors of N are

repeated, then there will be duplicate entries by using the above formula, so the number

of possible ways of representing N as sum of two squares will be lesser than that.

If N has certain number of ways to be represented as sum of two squares, then 2N can
also be represented as sum of two squares in the same number of ways. It becomes

9 9 9 9

evident from the following equation: 2 (a + b ) = (a + b) + (a - b) . Substituting with x
= 1 and y = 1 in the above equation, we simply arrive at this equation.

62

Sum of three squares: Every positive integer other than those of form 4 b (8c + 7), b > 0, c

>

can be represented as sum of three squares of no n- negative integers (including zero).

Sum of four squares: Every positive integer can be represented as sum of four squares of
non- negative integers (including zero). The Euler's identity clears this assertion.
(a 2 + b 2 + c 2 + d 2 ) (e 2 + f 2 + g 2 + h 2 ) = (ae + bf + eg + dh) 2 + (af - be + ch - dg) 2 + (ag -
ce + bh - df) 2 + (ah - de + bg - cf) 2

A.10 COMPUTING MINIMAL EQUAL SUMS OF LIKE POWERS

The Fermat's Last Theorem had been a long-standing unsolved problem in Mathematics
for over 300 years until it was resolved by Andrew Wiles in 1995, by using topics in very
advanced mathematics, elliptic curves and modular forms. It states that the Diophantine
equation x n + y n = z n has no solution in positive integers for any value of n > 3. Fermat
himself proved it for n = 4, by using method of infinite descent. Euler proved for n = 3,
and then Kummer proved it for all regular primes.

9 9 9

When n = 2, then the equation x + y = z becomes the Pythagoreas Equation, which has

9 9 9 9

infinitely many solutions in lowest terms given by ((m - n )/2, 2mn, (m + n )/2), where
m, n are odd, with m > n, gcd (m, n) = 1. This is the general solution for the above
Pythagoreas Equation for that (x, y, z) triplet in the lowest terms.

Beal's Conjecture: Analogous to the Fermat's Last Theorem, it has been conjectured that

the Diophantine equation x a + y b = z c has no solution in positive integers for (x, y, z)

when all of (a, b, c) > 3, when gcd (x, y, z) = 1.

If (x, y, z) are not co-prime to each other, then there exist solutions to the above equation

suchas756 3 + 945 3 =189 4 .

There also exist some solutions to the above equation when not all of (a, b, c) > 3, like

17 3 + 2 7 = 71 2 , 3 5 + 10 2 = 7 3 , ll 2 + 2 2 = 5 3 , 7 2 + 2 5 = 3 4 , 3 5 + ll 4 = 122 2 , 13 2 + 7 3 = 2 9 ,

n o ■)

1 +3 = 2 , etc. All these equations involve atleast some square term within them, for
example, 33 8 + 1549034 2 = 15613 3 , 1414 3 + 2213459 2 = 65 7 , 9262 3 + 15312283 2 = 113 7 ,
17 7 + 7627 1 3 = 21063928 2 , 43 8 + 96222 3 = 30042907 2 .

63

Computing Minimal Equal Sums of Like Powers: Let us now discuss about writing sum

of M k th powers as sum of N different k th power terms.

From the Fermat's Last Theorem, it follows that x 3 + y 3 = z 3 has no solution in positive

Q Q Q Q

integers. But for the equation x + y = z + w has many solutions.

10 3 + 9 3 = 12 3 + l 3 is Ramanujan's Taxicab Number. Also that 16 3 + 2 3 = 15 3 + 9 3 .

x + y + z = w also has got some solutions, too. For example, 3 +4 +5 = 6 .

That is for cubes. Coming to the fourth powers, x 4 + y 4 = z 4 has obviously got no solution
in positive integers. But that x 4 + y 4 = z 4 + w 4 as well as x 4 + y 4 + z 4 = w 4 both have got
atleast some solutions. 158 4 + 59 4 = 134 4 + 133 4 . 42248 1 4 = 414560 4 + 2175 19 4 + 95800 4 .

For the higher powers, I have got some examples with me as well.

144 5 = 133 5 +110 5 +84 5 +27 5

14132 5 +220 5 = 14068 5 +6237 5 +5027 5

23 6 +15 6 +10 6 = 22 6 +19 6 +3 6

966 8 +539 8 +81 8 = 954 8 +725 8 +481 8 +310 8 +158 8

3113 8 +2012 8 +1953 8 +861 8 =2823 8 +2767 8 +2557 8 +1128 8

Thus, it is conjectured that for forming such an equation of degree k, we need atleast k
terms.

A.ll SOME SPECIAL DIOPHANTINE EQUATIONS

A Diophantine Equation always seeks up for integer values as solutions, or that lattice

points within the Cartesian co-ordinate system.

The Diophantine equation 2 n - 7 = m 2 has got no solutions other than that for (m, n) =

(1,3), (3, 4), (5, 5), (11, 7), (181, 15).

This is closely related with Triangular Numbers of form 2 n - 1, since they are given by

the equation 2 n - 1 = m (m-l)/2, or that 2 n+3 - 7 = (2m - l) 2 .

Also, notice that whether the equation n! + 1 = m has got any other solutions besides that
for(m,n) = (5,4),(ll,5),(71,7).

64

The linear Diophantine equation is given by ax+by = c, and that equation has got
solutions for x, y if and only if c is a multiple of gcd(a, b).

Next, take up with that quadratic Diophantine equation, that is the Pell's Equation.

9 9

This equation Ix - 2y I = 1, has got solutions for x, y that are given by using the
convergents of the continued fraction expansion for \2. (x, y) is a solution for that above
equation if and only if x/y is a convergent of the continued fraction expansion of v2.
For example, the first few solutions are given as: (1, 1), (3, 2), (7, 5), (17, 12), (41, 29),
(99, 70)... Note that x; = 2x;.i + x;.2, and then that yi = 2yi_i + yi_2 as well. Also note that

9 9

the value of x - 2y alternates between 1 and -1 for each of these adjacent solutions.
In case it is not clear, 1/1, 3/2, 7/5, 17/12, 41/29, 99/70 are all convergents of continued
fraction expansion of \2, because of the fact that it is an approximation of value (as a
fraction) for this series, after every fixed number of terms within that:
1+1/(2+1/(2+1/(2+1/...)))

A.12 STRAIGHTEDGE AND COMPASS CONSTRUCTIONS

A straightedge is a ruler with no markings on it. Given a unit length, what lengths relative
to it can you construct by using straightedge and compass? What regular polygons can
you construct by using the given straightedge and compass? Let us answer that question
over here.

Given a straightedge and compass and a piece of paper, you are allowed to

Mark a random point on paper.

Draw a line, line segment or ray between two fixed points.

Draw a circle with center as one point, and passing through another existing point.

Mark intersection of two circles, two lines or line and circle.

Hide any unwanted points, lines or circles.

65

Marking middle of any two given points, drawing perpendicular bisector of any given

line segment, bisection of any given angle are possible with these facilities.

We can construct lengths of square roots of any integer. For example, \2 is the length of

hypotenuse of right angled triangle whose other two sides are 1 and 1, \3 is the length of

hypotenuse of right angled triangle whose other two sides are \2 and 1, and so on.

We can divide a line segment into any number of equal parts. Suppose that you want to

divide it into N equal parts, then take another parallel line segment containing N+l points

that are equally spaced, and then project that original line segment with this parallel line,

with the end points of the original line segment corresponding to that of the end points of

the parallel line segment.

Regular polygons: It is quite trivial to construct a polygon of 2N sides from a polygon of

N sides by bisecting the angle between the center of polygon to its adjacent vertices or by

joining the center of polygon to the mid-points of each of the edges of the polygon of N

sides and the points where these lines meet the circumscribed circle of the polygon of N

sides correspond to the other vertices of polygon of 2N sides.

To construct any regular polygon of N sides, you know the center of polygon and one of

its vertices, and you are required to construct the other vertices and then edges of the

polygon.

An equilateral triangle can be easily constructed by using straightedge and compass by

using the fact that its altitude is \3/2 times its side. Note that the altitude also bisects one

of its sides.

Construction of a square needs no explanation at all. It is very obvious enough.

A regular pentagon is not hard enough to construct. Every internal angle is 108°, and we

have that cos 72° = (a/5 - l)/4. Every angle between the center of pentagon to its adjacent

vertices is 72°, and we also have that the diagonal of the pentagon is pi = V(5-\5)/2 times

the radius of its circumscribed circle. In order to realize this, first of all we need to

construct a length of p2 = (v5 - l)/2, and then construct the length of pi = V(p2 2 + 1), by

using the Pythagoreas Theorem.

66

Regular hexagon can be constructed from equilateral triangle, regular octagon from

square, and so on.

A regular polygon of 15 sides can be constructed by using an equilateral triangle and then

a regular pentagon. Place all their vertices over a circle, with atleast one vertex of

equilateral triangle and regular pentagon being coincident. Then, there exist atleast one

vertex of the equilateral triangle and one vertex of regular pentagon that forms the

adjacent vertices of the regular polygon of 15 sides.

This trick cannot be applied to construct a regular polygon of 18 sides from an equilateral

triangle and regular hexagon, because of the fact that 3 and 6 are not co -prime.

In 1796, at the age of 19, Gauss realized that it is possible to construct a regular polygon
of 17 sides by using only straightedge and compass. A regular polygon of N sides can be
constructed by using straightedge and compass if and only if N is a power of 2 times 4 or
product of one or more distinct Fermat primes. The only known Fermat primes are 3, 5,
17, 257 and 65537. Then, it is known that there exist no more Fermat primes at all till F32.

A length X relative to a given unit length can be constructed by using straightedge and
compass if and only if X is an element of quadratic extension field of rational numbers. It
means that X should be a root of an irreducible polynomial equation with rational
coefficients whose degree is a power of 2. For example, cos {nil) is a root of irreducible
cubic 8x 3 +4x 2 -4x-l = 0, and that cos {nl9) is a root of irreducible cubic 8x 3 -6x+l = 0.
These lengths cannot be constructed because the degree of irreducible polynomial is 3,
and not a power of 2. Thus, it becomes evident that regular heptagon or regular nonagon
cannot be constructed by using straightedge and compass only. Note that there also exist
certain irrational numbers that do not satisfy any polynomial equation with rational
coefficients at all. These are called transcendental numbers, n and e are examples of
transcendental numbers. Other irrational numbers are said to be algebraic. Lengths of
transcendental numbers cannot be also constructed by using straightedge and compass
only.

67

We can construct a rational length times any given length, because you can make any
number of copies of given length, concatenate it, and then divide that line segment into
any number of equal parts.

Square root of any given length X is also possible. Take a line segment of length X+l,
and then draw a semicircular arc with this line segment as diameter. Then, draw a
perpendicular at a distance of 1 from the end point of its line segment. Then, the height at
which this perpendicular cuts the semicircular arc will be equal to "VX. In this way, we
can construct lengths like 4 \2 by using only straightedge and compass, that is the square
root of a/2 itself. (Figure to the left)

Suppose that you want to construct the lengths that are roots of quadratic equation x - ax
+ b = 0, then draw a circle with coordinates (a, b) and then (0, 1) as the end points of
diameter of the circle. Then, the X-coordinates of the points where this circle cuts up with
the X-axis will give away with the lengths that are the roots of this quadratic equation.
(Figure to the right)

(0,1)

I a (** °) -""" (*» °)

If we take up with a line segment with two parts of lengths 1 and b, and then project the
first part with length a, then the length of projection of second part will be ab
(Multiplication of two lengths)

If we take up with a line segment with two parts of lengths b and 1, and then project the
first part with length a, then the length of projection of second part will be a/b (Division
of two lengths) - (Illustrations being presented upon the next page).

68

Constructions that cannot be done with the aid of straightedge and compass only are:

Squaring a circle - Constructing a square with the equal area of any given circle. Here we

assume that circle of radius unit length. This involves construction of side of square with

length \7i. It is not possible since yn is a transcendental number.

Doubling a cube - Given a cube of unit length, construct a cube whose volume is twice

that of the original cube. This is not possible, since the length of edge of this new cube is

3 \2, which is a root of irreducible cubic x 3 - 2 = 0, with rational coefficients is of degree

3, that is not a power of 2 at all.

Trisection of an arbitrary angle is not possible by using only straightedge and compass as

well.

Some constructions that are possible by using only straightedge and compass include:
Tangent to circle from a point on circle, Tangents to circle from a point outside circle,
Trisection of right angle, A line that is parallel to any given line, perpendicular to a line
from a point on the given line, perpendicular to a line from a point outside the given line,
Angle bisectors, partitioning a given line segment into any number of equal parts,
centroid, orthocentre, circumcentre, circumcircle, incentre, incircle, excentre, excircle of
any given triangle, center of any given circle, diameter of circle, middle of two points,
perpendicular bisector of any given line segment, parallelogram, square with any given
side, square with any given diagonal, Translation of any given length on a line to any
other part over the same given line itself.

Some constructions that I am not sure of: Direct and inverse common tangents to any
given two circles, Ten cases of Appolonius problem, that is construction of circles that
are being tangential to any three given elements, where that given element may be either
a point, a line or an other circle itself, as well.

69

APPENDIX B
NUMERICAL RESULTS

B.l AURIFEUILLIAN FACTORIZATIONS

(2 4k " 2 +l) = (2 2k " 1 -2 k +l) (2 2k " 1 +2 k +l)

(3 6k " 3 +l) = (3 2k - x +l) (3 2k - x -3 k +l) (3 2k - x +3 k +l)

(5 10k - 5 -l) = (5^-1) (5 4k - 2 -5 3k - 1 +3.5 2k - 1 -5 k +l) (5 4k - 2 +5 3k - l +3.5 2k - l +5 k +l)

(6 12k " 6 +l) = (6 4k - 2 +l) (6 4k - 2 -6 3k - 1 +3.6 2k - 1 -6 k +l) (6 4k - 2 +6 3k - l +3.6 2k - l +6 k +l)

(7 14k - 7 +l) = (7 2k - x +l) a 6k - 3 -7 5k - 2 +3J 4k - 2 -7 3k - l +3J 2k - l -7 k +l) (l 6k - 3 +l 5k - 2 +3J 4k - 2 +l 3k - l +3J 2k -

J +7 k +l)

(10 20k-io +1) = (10 4k-2 +1) ( io^. l0 ^ +5A0 ^. 2 A0 5k - 2 +lA0 4k - 2 -2A0 3k ' l +5A0 2k - l AO k +l)

(10 8k - 4 +10 7k - 3 +5.10 6k " 3 +2.10 5k - 2 +7.10 4k - 2 +2.10 3k - 1 +5.10 2k - 1 +10 k +l)

(ll 22k - n +l) = (ll 2kl +l) (1 l 10k " 5 -l l 9k4 +5.1 l 8k4 -l l 7k3 -l l 6k3 +l l 5k2 -l l 4k2 -l l 3kl +5. 1 l 2kl -

U k +1) (11 10 k -5 +11 9k-4 +5 11 8 k -4 +11 7 k -3_ 11 6 k -3_ 11 5 k -2_ 11 4 k -2 +11 3k4 +5 n 2 k -l +llk+1)

(12 6k " 3 +l) = (12 2kl +l) (12 2k " 1 -2 2k - 1 .3 k +l) (12 2k " 1 +2 2k - 1 .3 k +l)

B.2 THE FIRST 33 FERMAT NUMBERS F n = 2 2 "+l

F = 3

Fi = 5

F 2 = 17

F 3 = 257

F 4 = 65537

F 5 = 641 x 6700417

F 6 = 274177 x 67280421310721

F 7 = 59649589127497217 x 5704689200685129054721

F 8 = 1238926361552897 x

93461639715357977769163558199606896584051237541638188580280321

F 9 = 2424833 x 7455602825647884208337395736200454918783366342657 x

741640062627530801524787141901937474059940781097519023905821316144415759

5047050080928 18711 693940737

70

Fio = 45592577 x 6487031809 x 4659775785220018543264560743076778192897 x

P 2 52

Fn = 319489 x 974849 x 167988556341760475137 x 3560841906445833920513 x P 564

F12 = 114689 x 26017793 x 63766529 x 190274191361 x 1256132134125569 x

568630647535356955169033410940867804839360742060818433 x C1133

F13 = 2710954639361 x 2663848877152141313 x 3603109844542291969 x

319546020820551643220672513 x C2391

F14 = 116928085873074369829035993834596371340386703423373313 x C 48 8o

F15 = 1214251009 x 2327042503868417 x 168768817029516972383024127016961 x

C98O8

Fie = 825753601 x 188981757975021318420037633 x C l9m

Fn = 31065037602817 x C 39 444

Fig = 13631489 x 81274690703860512587777 x C 78 884

Fi 9 = 70525124609 x 646730219521 x 37590055514133754286524446080499713 x

C157770

F20 = C3 15653

F 2 i = 4485296422913 x C 63 i294

F 22 = 64658705994591851009055774868504577 x C1262577

F 23 = 167772161 XC2525215

F24 = C5050446

F 25 has factors 25991531462657, 204393464266227713, 2170072644496392193

F 26 has a factor 76861124116481

F 27 has factors 151413703311361, 231292694251438081

F 28 has a factor 1766730974551267606529

F 29 has a factor 2405286912458753

F 30 has factors 640126220763137, 1095981164658689

F31 has a factor 46931635677864055013377

F 32 has a factor 25409026523137

71

B.3 FACTORS OF MERSENNE NUMBERS
(2 n -l, with n < 150 and n is prime)

n Prime Factorization of 2 n - 1

2 3

3 7
5 31
7 127
11 23.89
13 8191
17 131071
19 524287
23 47.178481

29 233.1103.2089

31 2147483647

37 223.616318177

41 13367.164511353

43 431.9719.2099863

47 2351.4513.13264529

53 6361.69431.20394401

59 179951.3203431780337

61 2305843009213693951

67 193707721.761838257287

71 228479.48544121.212885833

73 439.2298041.9361973132609

79 2687.202029703.1113491139767

83 167.57912614113275649087721

89 618970019642690137449562111

97 1 1447. 13842607235828485645766393

101 7432339208719.341117531003194129

103 2550183799.3976656429941438590393

72

107 162259276829213363391578010288127

109 745988807.870035986098720987332873

113 3391.23279.65993.1868569.1066818132868207

127 170141183460469231731687303715884105727

131 263.10350794431055162386718619237468234569

137 32032215596496435569.5439042183600204290159

139 5625767248687.123876132205208335762278423601

149 86656268566282183151.8235109336690846723986161

B.4 LIST OF MERSENNE, WAGSTAFF, REPUNIT PRIMES
Mersenne Primes are values of n such that 2" - 1 is prime

2 3 5 7 13 17 19 31 61 89 107 127 521 607 1279 2203 2281 3217 4253 4423
9689 9941 11213 19937 21701 23209 44497 86243 110503 132049 216091
756839 859433 1257787 1398269 2976221 3021377 6972593 13466917
20996011 24036583 25964951 30402457 32582657 37156667 42643801 43112609
Wagstaff Primes are values of n such that (2 n + l)/3 is prime
2 3 5 7 11 13 17 19 23 31 43 61 79 101 127 167 191 199 313 347 701 1709 2617
3539 5807 10501 10691 11279 12391 14479 42737 83339 95369 117239 127031
138937 141079 267017 269987 374321 986191 4031399
Repunit Primes are values of n such that (10 n - l)/9 is prime
2 19 23 317 1031 49081 86453 109297 270343

73

```